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14.  ABSTRACT 

Based  on  first  principles  supercell  approach,  we  have  successfully  initiated  and  continued  our  study  of 
defects  in  AIN,  GaN,  diluted  nitrides  and  related  wide  gap  materials.  We  have  also  investigated  and 
implemented  an  approach  to  calculate  the  x-ray  absorption  near  edge  structures  (XANES)  by  using  first 
principles  pseudopotential  in  conjunction  with  ab  initio  multiple  scattering  codes.  We  have  surpassed  our 
original  milestones  in  several  aspects  and  already  accomplished  important  results,  a)  We  have  calculated 
native  defects  in  AIN,  namely  VA1,  VN,  Ali,  and  Ni.  Most  of  the  results  are  in  agreement  with  previous 
work  by  Stampfl  and  Van  de  Walle.  Except  for  the  Ni,  where  we  have  shown  that  Ni  prefers  the 
split-interstitial  configuration,  b)  Potential  p-type  dopants  in  AIN,  Mg  and  Be,  have  been  studied.  We 
found  that  both  dopants,  when  substituting  for  Al,  give  a  rather  deep  level.  This  means  both  Mg  and  Be  are 
not  suitable  p-type  dopants  in  AIN.  c)  We  have  calculated  the  Ga  Frenkel  pairs  (interstitial  Ga  and  gallium 
vacancy  complexes)  in  GaN.  We  studied  both  the  stability  of  the  pair  at  different  separations  and  the 
barriers  for  the  pair  to  form/disintegrate.  Our  results  show  that  the  pair  at  ~  4  ?separation  is  stable  and 
has  a  binding  energy  of  about  3.3  eV.  However,  the  barrier  for  the  pair  to  recombine  is  only  ~  0.2  eV, 
making  it  unlikely  to  be  stable  at  room  temperature,  d)  In  dilute  nitrides,  the  nitrogen  atoms  are  highly 
strained.  Our  calculations  show  that  they  prefer  to  bind  with  small  impurity  atoms  such  as  C  or  Si  in  the 
form  of  CN  and  SiN  molecule,  substituting  on  the  As  site.  The  CNAs  gives  a  vibration  frequency  in  good 
agreement  with  the  local  vibration  mode  observed  in  an  infrared  spectroscopy  measurement.  The  SiNAs 
has  low  formation  energy  and  could  explain  the  mutual  passivation  of  shallow  donor  Si  and  isovalent  N  in 
dilute  nitride,  e)  We  investigated  the  local  structure  of  indium  oxynitride  at  various  0:N  ratios  using 
combined  experimental  x-ray  absorption  spectroscopy  and  first  principles  calculations.  Our  results  showed 
that  O  can  substitute  on  the  N  site  of  wurtzite  InN  upto  about  40%  then  the  alloy  started  to  have  mixed 
wurtzite  and  bixbyite  (In203  natural  structure)  phases.  The  simulated  x-ray  absorption  spectra  based  on 
the  proposed  models  are  in  a  very  good  agreement  with  the  measurements,  f)  We  have  revisited  the 
problem  of  the  H  binding  site  in  ZnO.  In  some  commercial  samples,  the  vibration  frequency  of  the  H  atoms 
indicates  that  they  are  binding  with  O  at  the  bond  center  site  whereas  in  some  samples  the  frequency 
indicates  that  they  are  at  the  anti-bonding  site.  Based  on  first  principles  calculations,  we  showed  that  Ca, 
found  only  in  some  samples,  plays  a  pivotal  role  in  switching  the  H  from  the  bond-center  site  to  the 
anti-bonding  site. 
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Executive  summary 

Based  on  first  principles  supercell  approach,  we  have  successfully  initiated  and 
continued  our  study  of  defects  in  AIN,  GaN,  diluted  nitrides  and  related  wide  gap 
materials.  We  have  also  investigated  and  implemented  an  approach  to  calculate  the  x-ray 
absorption  near  edge  structures  (XANES)  by  using  first  principles  pseudopotential  in 
conjunction  with  ab  initio  multiple  scattering  codes.  We  have  surpassed  our  original 
milestones  in  several  aspects  and  already  accomplished  important  results,  a)  We  have 
calculated  native  defects  in  AIN,  namely  Fai,  Fn,  Af  ,  and  N(.  Most  of  the  results  are  in 
agreement  with  previous  work  by  Stampfl  and  Van  de  Walle.  Except  for  the  N,,  where 
we  have  shown  that  N,  prefers  the  split-interstitial  configuration,  b)  Potential  p- type 
dopants  in  AIN,  Mg  and  Be,  have  been  studied.  We  found  that  both  dopants,  when 
substituting  for  Al,  give  a  rather  deep  level.  This  means  both  Mg  and  Be  are  not  suitable 
p- type  dopants  in  AIN.  c)  We  have  calculated  the  Ga  Frenkel  pairs  (interstitial  Ga  and 
gallium  vacancy  complexes)  in  GaN.  We  studied  both  the  stability  of  the  pair  at  different 
separations  and  the  barriers  for  the  pair  to  form/disintegrate.  Our  results  show  that  the 
pair  at  ~  4  A  separation  is  stable  and  has  a  binding  energy  of  about  3.3  eV.  However,  the 
barrier  for  the  pair  to  recombine  is  only  ~  0.2  eV,  making  it  unlikely  to  be  stable  at  room 
temperature,  d)  In  dilute  nitrides,  the  nitrogen  atoms  are  highly  strained.  Our 
calculations  show  that  they  prefer  to  bind  with  small  impurity  atoms  such  as  C  or  Si  in 
the  form  of  CN  and  SiN  molecule,  substituting  on  the  As  site.  The  CNAs  gives  a 
vibration  frequency  in  good  agreement  with  the  local  vibration  mode  observed  in  an 
infrared  spectroscopy  measurement.  The  SiNAs  has  low  formation  energy  and  could 
explain  the  mutual  passivation  of  shallow  donor  Si  and  isovalent  N  in  dilute  nitride,  e) 
We  investigated  the  local  structure  of  indium  oxynitride  at  various  0:N  ratios  using 
combined  experimental  x-ray  absorption  spectroscopy  and  first  principles  calculations. 
Our  results  showed  that  O  can  substitute  on  the  N  site  of  wurtzite  InN  upto  about  40% 
then  the  alloy  started  to  have  mixed  wurtzite  and  bixbyite  (In203  natural  structure) 
phases.  The  simulated  x-ray  absorption  spectra  based  on  the  proposed  models  are  in  a 
very  good  agreement  with  the  measurements,  f)  We  have  revisited  the  problem  of  the  H 
binding  site  in  ZnO.  In  some  commercial  samples,  the  vibration  frequency  of  the  H 
atoms  indicates  that  they  are  binding  with  O  at  the  bond  center  site  whereas  in  some 
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samples  the  frequency  indicates  that  they  are  at  the  anti-bonding  site.  Based  on  first 
principles  calculations,  we  showed  that  Ca,  found  only  in  some  samples,  plays  a  pivotal 
role  in  switching  the  H  from  the  bond-center  site  to  the  anti-bonding  site. 
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(Following  work  schedule  and  milestones  are  the  same  as  previously  submitted) 

Work  schedule: 


Month  0: 

o  Start  calculations  of  common  native  defects  and  dopants  (such  as  Mg  and  Be)  in 
AIN,  focusing  on  the  p- type  aspects  using  first  principles  ultrasoft  pseudopotential 
(USPP)  and  a  supercell  approach.  Since  both  Mg  and  Be  are  theoretically  reported 
to  be  good  p- type  dopants  in  GaN,  they  are  good  candidates  for  AIN  and  GaAIN 
alloys. 

o  Calculations  of  some  impurities  in  diluted  nitrides  (such  as  C-N  complexes  and  Si- 
N  complexes  in  GaAsN)  using  first  principles  supercell  approach. 

Month  6: 

o  Low  energy  native  defects  in  AIN  and  the  potential  dopants  for  p- type  AIN  will  be 
further  investigated.  Their  properties  in  AlGaN  alloys  will  be  studied. 

o  Start  calculations  of  Frenkel  pairs  in  GaN  by  first  principles  supercell  approach. 
The  stability  of  the  pairs  will  be  investigated  by  varying  the  pair  distance.  If 
needed,  effects  of  band  gap  corrections  might  be  investigated  by  using  LDA+t/ 
approach. 

o  Simulating  the  signatures  of  interesting  defects  in  comparison  with  available 
measurements  (such  as  XAS,  XPS,  PL,  DLTS,  EPR,  or  IR  spectroscopy). 
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Milestones: 


6  month: 

o  Obtained  the  results  of  at  least  four  native  defects  (and/or  impurities)  in  AIN. 

o  Obtained  some  stable  impurity  structures  and  their  properties  (such  as  vibration 
signatures)  in  diluted  nitrides.  For  examples,  CN  complexes  in  GaAsN  and  SiN 
complexes  in  GaAsN. 

12  month: 

o  Complete  results  on  simple  native  defects  as  well  as  Mg-  and  Be-doped  AIN. 

o  Obtained  the  binding  energies  and  other  properties  of  the  Ga,  -  Foa  and  N,  -  Fn 
Frenkel  pairs  in  GaN. 

o  Results  on  the  binding  energies  and  stabilities  of  Frenkel  pairs  as  a  function  of  the 
pair-distance. 

o  Signatures  (that  can  be  compared  with  experimental  measurements)  of  some 
interesting  defects  and  dopants  in  AIN  (and  AlGaN  alloys). 

o  Implications  of  some  interesting  AIN  defects  when  formed  in  AlGaN  alloys. 
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Research  Activities 

Computational  codes 

First-principles  calculation  approach  based  on  density-functional  theory  and  ab  initio 
ultrasoft  pseudopotentials  (USPP)  is  the  main  approach  for  this  work.  We  used  the  codes 
called  Vienna  Ab-initio  Simulation  Package  (VASP)  [1],  For  the  x-ray  absorption 
spectroscopy  (XAS)  simulation,  we  use  the  ab  initio  multiple  scattering  code  called  FEFF 
(version  8.2)  [2],  Note  that  before  calculating  the  XAS,  the  accurate  atomic  relaxations 
are  performed  using  the  first  principles  USPP  codes  (VASP).  The  detail  of  the  parameters 
and  approximations  used  in  the  calculations  is  slightly  varied,  depending  on  specific 
material.  Further  specific  computational  detail  for  each  system  under  study  can  be  found 
in  each  attached  publications. 

A  brief  description  of  the  defect  theory 

Calculations  for  defects  in  semiconductor  are  calculated  using  a  supercell 
containing  32  to  96  atoms  and  a  set  of  non-T  k-points  sampling. 

In  supercell  calculations,  the  formation  energy  of  a  defect  X  in  charge  state  q  is 

£/[*']  =  E„x(X- ) -  Emx(bu!k)  +  n„K  +  qEP,  [  1  ] 

where  Etot.sdX1)  is  the  total  energy  of  a  supercell  containing  the  defect,  EmSc(bulk)  is  the 
total  energy  of  the  same  supercell  but  without  any  defect.  na  is  the  number  of  the  atoms 
of  specie  a  being  removed  from  (positive  value)  or  added  to  (negative  value)  a  bulk 
supercell  to  form  defect  X.  /j,a  is  the  reservoir  energy  of  specie  a  (the  chemical 
potential).  £>  is  the  Fermi  level,  referenced  to  the  valence-band  maximum  of  bulk,  Ev. 

The  chemical  potentials  in  Eq.  [1]  depend  on  the  growth  conditions.  For  example, 
(in  the  case  of  defects  in  AIN)  under  Al-rich  case,  jUm  =  formation  energy  of  solid  metal 
Al.  In  order  for  AIN  to  be  stable  during  the  growth,  we  also  require  that  //Ai  +  jUn  =  //Ain 
=  constant  <F>.  This  fixes  //n  =  ®  -  //ai-  The  chemical  potential  of  other  elements  are 
referenced  to  their  natural  phase.  For  example  the  chemical  potential  of  C,  jUc  is 
referenced  to  the  diamond  limit,  i.e.  /Jc  =  £’t0t(diamond)/2.  The  total  energy  of  diamond  is 
divided  by  two  because  there  are  two  carbon  atoms  per  a  diamond  unit  cell. 
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Native  defects  in  AIN 


We  followed  the  computational  approach  that  we  previously  used  to  study  defects 
in  GaN  [3]  and  Ti02  [4]  to  study  native  defects  in  AIN.  Test  calculations  show  that  using 
the  supercell  size  of  36  or  96  atoms  result  in  similar  formation  energy;  indicating  that  96- 
atom  supercell  is  sufficiently  large  for  the  defect  calculations  of  wurtzite  AIN.  This  is 
similar  to  the  case  of  wurtzite  GaN.  Native  defects  in  AIN  have  been  previously  studied 
by  Stampfl  and  Van  de  Walle  [5]  using  the  pseudopotential  planewave  method.  Their 
calculations  use  mainly  a  32-atom  zincblende  supercell  except  some  calculations  using 
the  wurtzite  supercell  or  96-atom  wurtzite  supercell  for  convergence  tests.  Here  we  use 
slightly  different  computational  approach  (USPP  planewave);  allowing  the  calculations 
with  larger  wurtzite  supercell  (all  calculations  are  done  with  96-atom  cell).  Most  of  the 
results  are  in  a  reasonable  agreement  with  the  previous  work  by  Stampfl  and  Van  de 
Walle  [5]  with  some  exceptions  as  will  be  discussed  below. 

1)  Vacancies  (aluminum  and  nitrogen  vacancies,  Vai  and  Vn) 

Similar  to  many  semiconductors,  vacancies  are  the  low  energy  defects  in  AIN.  As  shown 
in  Fig.l,  under  Al-rich  conditions,  it  is  more  difficult  to  create  aluminum  vacancy  (Vai); 
leading  to  higher  energy  of  Vai.  Since,  under  that  condition,  nitrogen  is  scarce,  the 
nitrogen  vacancies  can  be  created  easily;  leading  to  low  V n  energy.  Under  N-rich 
conditions,  the  formation  energy  of  Vai  shifted  down  and  Vn  shifted  up.  For  Vai,  the 
possible  charge  states  varied  from  neutral  near  the  VBM  to  3-  at  higher  Fermi  energy 
above  1.5  eV  (The  charge  state  can  be  determined  from  the  slope  of  the  plot.)  As  a 
result,  Vai  can  be  considered  as  a  very  deep  acceptor  that  can  accept  up  to  three  electrons. 
For  Vn,  there  are  two  possible  charge  states,  the  3+  charge  state  at  the  Fermi  energies 
near  the  VBM  and  the  1+  charge  state  for  the  higher  Fermi  energies.  The  sharp  drop  in 
formation  energy  at  Fermi  energy  above  3.5  eV  is  an  artifact  in  the  calculations  from  the 
occupations  of  the  calculated  conduction  band  minimum  that  have  to  be  ignored.  This 
issue  has  also  been  discussed  in  Ref.  [5].  (In  the  calculations,  the  band  gap  appears  too 
small  due  to  the  well-known  DFT  gap  underestimation.)  In  reality  the  1+  charge  state 
would  extended  along  the  dotted  line  all  the  way  to  the  CBM.  As  a  result,  Vn  is  an  active 
donor  that  can  donate  three  electrons  when  the  Fermi  energy  is  near  the  VBM  (i.e.,  when 
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the  sample  is  p- type)  or  donate  one  electron  when  the  Fermi  energy  lies  in  the  range  from 
1  eV  to  the  bottom  of  CBM  (n-type  condition).  These  vacancy  defects  could  play 
important  roles  in  compensating  the  dopants.  For  p- type  doping,  Fn,  as  a  triple  donor, 
would  play  an  important  role  in  compensating  the  acceptor  dopants.  To  avoid  such 
compensation,  one  should  perform  p- type  doping  under  N-rich  conditions  where  the 
formation  energy  of  Fn  is  the  highest.  On  the  other  hand,  for  »-type  doping,  VA\,  as  a 
triple  acceptor,  would  play  an  important  role  in  compensating  the  donor  dopants. 
However,  Fai  is  only  active  under  n-type  conditions.  By  itself  F\i  does  not  act  as  an 
acceptor  because  the  s(0/-)  is  very  deep,  i.e.  the  defect  level  is  far  away  from  the  VBM. 
This  is  in  contrast  to  Fn  which  always  act  as  a  donor  regardless  of  the  Fermi  level  of  the 
sample 


Al-rich  N-rich 


Figure  1:  Calculated  formation  energy  as  a  function  of  the  Fermi  level  for  vacancy  defects 
under  the  Al-rich  (left)  and  N-rich  (right)  limits.  The  Fermi  level  is  references  to  the  top  of 
the  AIN  VBM. 
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2)  Interstitials  (interstitial  aluminum  and  nitrogen,  Alm,  and  Nu„) 

Interstitial  defects  are  generally  high  in  energy  than  vacancy  defects,  especially  for 
covalent  materials.  However,  in  ionic  materials  the  charged  interstitial  defects  can  have 
low  formation  energies  due  to  the  columbic  energy  gain  between  the  cation-anion 
bindings.  For  example,  interstitial  titanium  in  4+  charge  state  has  been  found  to  have  low 
formation  energy  in  TiCh  [4],  Here,  we  found  that  Al,„,  in  3+  charge  state  has  low 
formation  energy  in  AIN  under  Al-rich  conditions.  From  the  formation  energy  plot  in 
Fig.  2,  Al;>„  is  stable  in  3+  charge  state  for  all  Fermi  energy;  indicating  that  it  is 
exclusively  a  donor.  Note  that,  the  drop  in  formation  energy  near  CBM  is  an  artifact  in 
the  calculations  from  the  occupations  of  the  calculated  conduction  band  minimum  that 
have  to  be  ignored.  Interstitial  nitrogen  in  AIN  is  stable  in  the  split  interstitial  form 
which  can  be  viewed  as  a  nitrogen  molecule  substituting  the  N  site  [3,  6],  N„„  in  this 
form  is  lower  in  energy  than  an  isolated  form  by  at  least  2  eV.  Since  Stampfl  and  Van  de 
Walle  [5]  did  not  investigate  N;>„  in  this  form,  their  formation  energy  of  this  defect 
appears  higher.  Regarding  the  charge  state,  Nm/  is  amphoteric  defect,  i.e.  can  behave  as 
both  donor  and  acceptor  to  compensate  for  the  leading  carrier  in  the  sample.  However,  it 
has  reasonably  low  energy  only  under  N-rich  conditions  and  Fermi  energy  near  the  CBM. 
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Figure  2:  Calculated  formation  energy  as  a  function  of  the  Fermi  level  for  intersitial  defects 
under  the  Al-rich  (left)  and  N-rich  (right)  limits.  The  Fermi  level  is  references  to  the  top  of 
the  AIN  VBM. 
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3)  Antisites  (aluminum  and  nitrogen  antisites,  A1n  and  Naj) 

Generally  antisite  defects  have  high  formation  energy,  especially  in  ionic  materials.  For 
AIN,  which  is  rather  ionic,  the  antisite  defects  are  reported  to  be  higher  in  energy  than 
vacancy  defects  [5],  However,  AIN  has  rather  large  heat  of  formation.  Our  calculated 
value  of  AH  =  3.6  eV  is  consistent  with  Ref.  [5]  of  3.5  eY.  This  gives  a  wide  variation 
of  (theoretical)  possible  growth  conditions  from  Al-rich  to  N-rich.  For  vacancy  defects, 
this  means  the  variation  in  formation  energy  of  3.59  eV  from  Al-rich  to  N-rich. 
However,  for  antisite  defects,  this  means  the  variation  of  3.6  x  2  =  7.2  eV.  As  a  result, 
the  calculated  formation  energies  of  these  antisite  defects  suggested  that  under  extreme 
Al-rich  and  N-rich  conditions  they  can  be  stable.  As  shown  in  Fig.  3,  the  A1N  could  be 
stable  in  4+  charge  state  under  Al-rich  and  p- type  conditions  whereas  the  Nai  could  be 
stable  in  4-  charge  state  under  N-rich  and  n-type  conditions.  Our  calculations  show  that 
these  two  highly  charged  antisite  defects  can  be  stable  over  other  native  defects  for  a 
small  range  of  Fermi  energy  under  these  two  extreme  Al-rich  and  N-rich  conditions. 
This  is  slightly  different  from  Stampfl  and  Van  de  Walle  results  [5]  that  found  the  antisite 
defects  slightly  higher  in  energy  than  other  native  defects  and  have  never  been  the  lowest 
energy  defects  in  any  conditions.  This  is  due  partly  to  the  fact  that  our  calculated  heat  of 
formation  is  slightly  higher  than  theirs;  allowing  our  extreme  Al-rich  and  N-rich 

conditions  to  extend  a  little  bit  wider. 

Al-rich  N-rich 


Figure  3:  Calculated  formation  energy  as  a  function  of  the  Fermi  level  for  antisite  defects 
under  the  Al-rich  (left)  and  N-rich  (right)  limits.  The  Fermi  level  is  references  to  the  top  of 
the  AIN  VBM. 
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Magnesium  and  beryllium  doped  AIN 

Magnesium  has  been  used  successfully  to  dope  GaN  into  p- type;  leading  to 
commercial  blue  lasers  in  today  devices.  We  have  previously  suggested  that  beryllium  is 
a  potentially  good  p- type  dopant  in  GaN  as  well  [7],  However,  beryllium  can  potentially 
form  interstitial  (Be,„d  which  is  a  double  donor.  Moreover  the  interstitial  is  also  likely  to 
bind  with  the  substitutional  acceptor  (BeGa)  turning  the  acceptor  into  a  single  donor.  In 
our  previous  work  [7],  we  suggested  that  if  one  can  overcome  interstitial  Be  problem, 
beryllium  is  a  potentially  better  acceptor  than  Mg.  Because  the  lattice  constant  of  AIN  is 
very  similar  to  that  of  GaN,  both  Mg  and  Be  are  potential  candidates  as  the  p- type 
dopants  in  AIN.  However,  experiment  results  show  that  Mg  is  a  deep  acceptor  in  AIN 
(with  ionization  energy  of  around  0.4  -  0.6  eV)[8].  Here,  we  systematically  study  both 
Mg  and  Be  as  dopants  in  AIN.  Interstitial  Mg  and  Be  (Mg,,,,  and  Be,„f)  are  also  study  to 
access  the  likelihood  that  they  can  cause  spontaneous  compensation. 

1)  Magnesium  substituting  for  Al  and  magnesium  interstitial  (MgAi  and  Mg  jnt). 

Figure  4  shows  the  formation  energies  of  Mg  substituting  for  Al  and  interstitial 
Mg  in  AIN.  It  is  similar  to  the  case  of  Mg  in  GaN  that  the  Mg  interstitial  always  a  double 
donor,  i.e.  exist  only  in  charge  state  2+  regardless  of  Fermi  level.  Here  the  formation 
energy  of  Mg/„,  under  Al-rich  and  N-rich  are  different  because  the  chemical  potential  of 
Mg  is  calculated  differently  for  each  case.  (If  the  chemical  potential  of  magnesium  is 
calculated  from  Mg  metal  the  formation  energy  of  Mg„„  would  not  depend  on  whether 
the  condition  is  Al-rich  or  N-rich.)  For  Al-rich  conditions,  the  chemical  potential  of  Mg 
is  limited  to  the  formation  of  (fee)  Mg  metal.  However,  for  N-rich  conditions,  it  is  the 
Mg3N2  formation  that  precipitates  Mg  chemical  potential.  As  a  result,  the  formation 
energy  of  Mg,,,,  under  Al-rich  condition  is  slightly  lower  due  to  the  higher  possible  Mg 
chemical  potential.  For  the  substitution  Mg,  obviously  the  N-rich  condition  provides  the 
best  opportunity  for  Mg  to  substitute  for  Al.  Therefore,  the  formation  energy  of  MgAi  is 
lower  for  N-rich  condition.  According  to  our  calculations,  MgAi  is  a  rather  deep  acceptor 
with  the  ionization  energy  of  0.69  eV.  Under  Al-rich  condition,  Mg„„  has  rather  low 
formation  energy  and  can  compensate  the  MgAi  if  there  is  no  other  lower  energy  donor. 
This  suggests  that  it  is  difficult  to  dope  Mg  into  AIN  without  forming  an  equal  amount  of 
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compensating  Mg,-B(  under  Al-rich  conditions.  However,  one  might  success  in  doping 
Mg  under  N-rich  conditions  where  the  Mg;>„  has  higher  formation  energy.  Nevertheless, 
MgA 
temp 

Al-rich  N-rich 


Fermi  energy  [eV]  Fermi  energy  [eV] 


Figure  4:  Calculated  formation  energy  as  a  function  of  the  Fermi  level  for  MgAi  and  Mgjn/ 
under  the  Al-rich  (left)  and  N-rich  (right)  limits.  The  Fermi  level  is  references  to  the  top  of 
the  AIN  VBM. 


2)  Beryllium  substituting  for  Al  and  beryllium  interstitial  (BeAi  and  B  Qint). 

The  formation  energies  of  Be  substituting  for  Al  and  interstitial  Be  in  AIN  is 
similar  to  the  cases  of  Mg  in  AIN  shown  in  the  previous  paragraph.  The  Be  interstitial  is 
also  a  double  donor.  (The  neutral  charge  near  the  top  of  CBM  is  probably  an  artifact  in 
the  calculation  caused  by  the  underestimation  of  the  band  gap  which  leads  to  the  wrong 
occupation.)  Figure  5  shows  the  formation  energy  of  BeAi  and  Be,„,  in  the  same  fashion 
as  Fig.  4  which  is  the  case  of  Mg.  For  Al-rich  conditions,  the  chemical  potential  of  Be  is 
limited  to  the  formation  of  (fee)  Be  metal.  However,  for  N-rich  conditions,  it  is  the 
Be3N2  formation  that  precipitates  Be  chemical  potential.  In  comparison  with  the  case  of 
Mg  in  AIN,  the  formation  energy  of  interstitial  Be  is  much  lower  than  that  of  Mg;>„.  This 
suggests  that  Be„;/  can  form  much  easier  than  Mg,-,,,  which  would  cause  more 
compensating  problems  in  AIN.  Note  that,  Be  interstitial  is  also  predicted  to  cause  more 
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compensating  problem  than  Mg  in  the  case  of  GaN  [7].  The  N-rich  condition  provides 
the  best  opportunity  for  Be  to  substitute  for  Al.  Therefore,  the  formation  energy  of  BeAi 
is  lower  for  N-rich  condition.  Unfortunately,  according  to  our  calculations,  BeAi  is  also  a 
rather  deep  acceptor  with  the  ionization  energy  of  0.67  eV  (almost  the  same  as  MgAi). 
Under  Al-rich  condition,  Be,,,,  has  low  formation  energy  and  can  compensate  the  BeAi  if 
there  is  no  other  lower  energy  donor.  Therefore,  it  is  difficult  to  dope  Be  into  AIN  and 
avoid  compensating  BQint  under  Al-rich  conditions.  N-rich  is  a  better  condition  for 
doping  Be  mainly  because  Be,„,  has  higher  formation  energy.  However,  the  fact  that  BeAi 
has  high  ionization  energy  makes  it  a  poor  p- type  dopant  in  AIN. 


Figure  5:  Calculated  formation  energy  as  a  function  of  the  Fermi  level  for  BeAi  and  Be,,„ 
under  the  Al-rich  (left)  and  N-rich  (right)  limits.  The  Fermi  level  is  references  to  the  top  of 
the  AIN  VBM. 
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Interstitial  gallium  and  gallium  vacancy  complexes  (Frenkel  pairs)  in  GaN 

Native  defects  can  agglomerate  into  complexes  during  the  growth  process,  or 
subsequently  cool  down  after  post-annealing.[9-ll]  Their  existences  can  affect  the 
intentional  dopants  or  the  electrical  property  to  the  samples. [12]  One  of  the  direct 
approaches  to  create  large  amount  of  native  defects  is  the  electron  irradiation  process. 
After  electron  irradiation  process,  interstitials  and  the  corresponding  vacancies  are 
expected.  Depending  on  the  energy  of  the  irradiated  electron,  Ga  or  N  can  be  selectively 
knocked  out  of  their  sits.  There  are  many  experimental  works  on  the  electronic  and 
optical  properties  of  samples  after  electron  irradiation  in  GaN.  The  attempts  to  detect 
these  defects  have  been  carried  out. [13- 17]  However,  theoretical  studies  on  these 
complexes  are  scarce.  In  this  program,  we  focused  our  attentions  on  a  complex  which 
consists  of  an  interstitial  gallium  and  gallium  vacancy  ( Ga/  -  VGa ),  so-called  “gallium- 

Frenkel  pair”.  These  defects  are  expected  to  form  under  high  energy  electron  irradiation. 
Here,  we  have  performed  first  principles  supercell  calculations  of  gallium-Frenkel 
complex  as  well  as  relevant  defects  in  wurtzite  GaN.  The  possible  configurations  and  the 
migration  mechanism  as  well  as  the  corresponding  energy  barriers  are  investigated.  For 
the  migration  study,  the  nudge  elastic  band  (NEB)  method  [18]  was  used. 

Configurations  of  gallium  Frenkel  pairs  ( Ga(  -VGa) 

The  parent  defects  of  gallium-Frenkel  pair  are  interstitial  Ga  (Ga,)  and  Ga 
vacancy  (Toa)-  Gallium  interstitial  can  occur  in  3+,  2+,  and  1+  charge  states,  depending 
on  the  Fermi  energy  of  the  sample.  In  its  3+  charge  state,  Ga:’+  favors  a  site  near  the 

octahedral  interstitial  site.  This  site  locates  along  the  center  of  hexagonal  channel  which 
is  called  O'  site  [19],  This  O'  is  only  0.16  A  away  from  an  ideal  O  site  (the  symmetric 
site  by  symmetry).  At  O'  site  the  distance  to  the  six  Ga  nearest  neighbors  are  roughly 
the  same  (to  within  0.15  A).  The  defect  level  of  Ga,  is  found  to  lie  in  the  upper  part  of 
the  band  gap  and  has  ^-character.  The  level  is  completely  occupied  in  1+  charge  state 
and  empty  in  3+  charge  state.  Gallium  vacancy  can  occur  in  3-,  2-,  1-,  and  neutral 
charge  states.  The  defect  levels  are  raised  up  from  the  nitrogen  dangling  bonds, 
combined  into  ^-states.  The  ^-states  are  split  into  a  doublet  and  singlet  states  [20]  where 
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the  singlet  states  locate  under  valence-band  maximum  (YBM),  thus  are  already  occupied 
and  the  doublet  state  lie  near  above  VBM.  In  the  neutral  charge  state  the  doublet  states 
have  only  one  electron  occupied.  They  can  be  filled  with  up  to  three  more  electrons. 


Figure  6:  Atomic  structure  of  the  lowest  energy  Ga-Frenkel  pair  (1st  Garl'(i;1)  from  our 
calculations:  (a)  before  relaxation  and  (b)  after  relaxation.  Top  panels  show  a  top  view  and 
lower  panels  a  perspective  view. 


Because  Ga,  and  Vq.a  have  opposite  charge,  they  have  Coulom  attraction  to  each 
other.  However,  if  the  two  defects  are  placed  next  to  each  other,  they  are  spontaneously 
recombined  (Ga,  moves  into  the  center  of  Vqa  and  the  two  annihilated).  Even  if  Ga,  is 
placed  at  the  next  nearest  neighbor  of  Vg&,  Ga,  can  still  spontaneously  diffuse  via  the 
kick-out  mechanism  and  the  defects  are  annihilated.  The  first  stable 
Ga,.  -  fGa  configuration  occurs  when  the  pair  is  separated  by  ~  4  A  along  the  hexagonal 

channel  in  which  Ga,  is  placed  at  octahedral  site  as  shown  in  Fig.  6  (a).  Our  calculations 
show  that  all  attempts  to  place  with  two  defects  with  the  separating  distance  of  less  than 
this  would  lead  to  spontaneous  annihilation.  The  configuration  in  Fig.  6  is  labeled 
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“1st  Ga;  -FGa  ”,  where  1st  is  denoted  for  the  lowest  formation  energy.  We  found  that  after 

relaxation  the  Ga,  has  moved  by  0.56  A  along  the  [-1100]  direction  and  subsequently 
impels  the  two  nearest  Ga  atoms  toward  the  vacancy  with  a  large  displacement  of  1.17  A 
(see  Fig.  6(b)).  In  the  process,  a  nitrogen  atom  which  has  one  broken  bond,  resulting  from 
removal  of  a  Ga  atom,  relaxed  outward  significantly. 


Figure  7:  Atomic  structure  of  the  Ga-Frenkel  pair  in  2st  Ga,- FCa  configuration  from  our 
calculations:  (a)  before  relaxation  and  (b)  after  relaxation.  Top  panels  show  a  top  view  and 
lower  panels  a  perspective  view. 


The  next  stable  configuration  with  the  higher  energy,  labeled  2nd  Ga;  -  VGa  is 

found  when  the  Ga,  is  placed  at  the  next  nearest  octahedral  site  along  the  hexagonal 
channel  with  a  distance  of  6.4  A  away  from  the  ideal  Voa  position,  as  shown  in  Fig  7  (a). 
After  relaxation,  the  Ga,  spontaneously  move  toward  three  nitrogen  atoms  with  a  distance 
of  0.88  A  along  [-1100]  direction,  result  in  repulsion  to  the  nearby  Ga  atom  which  is 
moved  off  its  site  by  ~  1 .3  A  toward  the  closest  nitrogen  atom  with  a  broken  bond. 
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Beyond  the  2nd  Ga(  -  VGa ,  we  have  also  searching  for  the  next  non-recombined 
configurations.  However,  we  found  that  the  3rd  and  4th  non-recombined  Ga;  -  VGa  pairs 
have  almost  no  binding  as  will  be  illustrated  next. 


Energies  of  gallium  Frenkel  complexes. 

Figure  8  show  the  calculated  formation  energies  as  a  function  of  Fermi  energy  in 
Ga-rich  conditions.  Note  that  the  total  formation  energy  of  Ga;  -  F0a  pair  would  be  the 
same  for  both  Ga-rich  and  N-rich  conditions. 


Figure  8:  Formation  energies  as  a  function  of  Fermi  energy  for  Ga„  FG;l,and  Ga-Frenkel  pair 
(Gar  ('i,.,)  under  Ga-rich  condition.  Fermi  energy  is  referenced  to  the  AIN  VBM.  The  labels 
n,h  where  n  =1,  2,  3,  4  are  referred  to  the  lowest  to  highest  formation  energy  of  the  Frenkel 
pairs  which  is  also  corresponding  to  the  pair  separation  from  smallest  to  the  largest  distance. 
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The  calculated  defect  levels  of  isolated  Ga,  and  FGa  are  similar  to  those  reported 
in  Ref  [19,  20].  Ga  interstitial  can  donate  up  to  3  electrons.  As  the  Fermi  energy  is  raised 
the  5-like  defect  state  near  the  CBM  can  be  populated  by  two  electrons;  turning  the  defect 
into  1+  charge  state.  Note  that,  the  2+  charge  state  is  never  stable  (Ga,  is  a  negative- U 
system).  Ga,  acts  exclusively  as  a  donor.  The  defect  level  of  Ga,  is  found  at 
approximately  2  eY  above  the  VBM;  leading  to  the  transition  level  s{ 3  +  /l+)«2  eV 
above  the  VBM. 

In  the  neutral  charge  state  FGa  is  triply  occupied,  with  defect  levels  lie  close  to  the 
VBM.  These  levels  are  derived  from  nitrogen  /^-states  which  have  t2  symmetry  and  can 
accept  three  more  electrons;  making  FGa  a  (deep)  triple  acceptor.  As  the  Fermi  energy  is 
raised,  the  defect  states  can  be  occupied  one-by-one  resulting  in  1-,  2-  and  3-  charge 
states  [20],  These  levels  lead  to  the  transition  levels  e(0/~),  s(-/  2-),  and  £•(2-/3-) 
of  0.2  eV,  0.5  eV  and  0.9  eV,  respectively. 

The  calculated  formation  energies  for  the  four  considered  configurations  of 
Ga;. -FGa  pair  as  a  function  of  Fermi  energy  are  shown  in  Fig.  8.  The  formation  of 

Ga,  -FGa  follows  the  reaction:  Ga(  +FGa  — » (Ga ,  -FGa)  where  the  binding  energy  (F’bind) 
is  defined  as  the  difference  between  the  sum  of  Ga,  and  FGa  energies  (shown  as  dashed 
lines  in  Fig.  8)  and  the  energy  of  Ga;  -  FGa .  For  neutral  Ga/  -  VGa ,  the  binding  energy  is 

defined  as  Ebind(G af  -  FGa)  =  Ef[(Gai  - FGa)°] - Ef[ Ga^+] - Ef[V^] .  We  find  that  only 
the  first  two  lowest  energy  configurations  (1st  and  2nd  Ga/  -VGa)  give  reasonably  large 

binding  energies  ~  3.3  and  2.4  eV.  The  other  two  configurations  (with  larger  separation) 
have  very  little  binding  energy  (~  0.1  eV);  implying  that  they  are  not  stable.  The  large 
binding  energy  of  the  1st  and  2nd  Ga(  -FGa  configurations  indicates  that  they  are  stable 

against  the  decomposition  into  isolate  parent  defects.  (However,  whether  they  are  stable 
against  annihilation  into  a  perfect  crystal  is  subject  to  further  investigations.)  The  large 
binding  energies  of  1st  and  2nd  Ga;  -  FGa  configurations  can  be  explained  by  Coulomb 

attraction  between  the  opposite  charges  of  Ga;’+  and  VGa . 
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Valence  band 

The  schematic  illustration  for  electronic  level  in  neutral  charge  state  induced  by 
the  1st  Ga;.  -  VGa  configuration  is  shown  in  Fig  9.  This  indicates  that  the  defect  complex 
is  electrically  inactive  at  almost  entire  Fermi  energy  range  except  near  the  band  edges. 
The  2nd  Ga;  -  V(a  configuration  exhibits  similar  electronic  levels  with  slight  variation  in 
the  positions  (not  shown). 

Formations  and  separations  of  gallium  Frenkel  complex. 

To  investigate  the  energy  of  Ga/  -  VGa  as  a  function  of  the  pair  separation,  we 

calculated  the  energy  profile  of  the  pair  separation  by  using  nudge  elastic  band  method 
(NEB).  The  calculation  is  started  with  a  defect-free  GaN  bulk  supercell  as  a  ground-state 
configuration.  Then,  a  Ga  atom  is  dislocated  out  from  its  site  to  the  nearest  interstitial 
site;  leaving  a  vacancy  behind.  Next,  the  Ga,  atom  is  migrate  further  to  form  the 
1st  Ga;  -VGa  configuration  by  the  kick-out  mechanism  [19].  The  energy  profile  for  the 
migration  is  mapped  by  using  the  NEB  method  [18].  As  the  next  step,  the  energy  profile 
in  the  further  migration  from  the  1st  Ga;  -  VGa  configuration  to  the  2nd  Ga,-  Tea 
configuration  is  mapped  out.  It  is  found  that  the  energy  needed  for  migration  from 
ground  state  to  the  1st  Ga,- Eoa  configuration  is  ~  8.33  eV  which  is  the  formation  energy 
of  1st  Ga,-F(ja  configuration.  Additional  energy  of  0.86  eV  is  needed  to  further  transform 
into  the  2nd  Ga,-Eoa  configuration.  Overall,  this  entire  profile  shows  that  the  1st  (as  well 
as  2nd)  Gaz-Eca  configuration  is  metastable  and,  in  principle,  can  be  stabilized.  However, 

they  can  be  recombined  if  the  Ga,  atom  gains  sufficient  energy  to  overcome  the  barriers 
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Figure  9:  Schematic  illustration  of  the 
defect  levels  induced  by  the  1st  Ga,-  V(:a 
configuration  in  the  GaN  band  gap.  The 
solid  dots  show  the  electron  occupation 
for  the  neutral  charge  state. 


which  are  0.17  eV  and  0.04  eV  for  the  1st  and  2nd  Gar  Foa  configuration,  respectively. 
Since  Ga,  is  known  to  be  very  diffusive  at  room  temperature  [3,  16],  the  Ga;-  Vc,a  is 
unlikely  to  be  stable  unless  the  sample  is  kept  at  low  temperature.  This  energy  profile 
could  be  used  for  further  study  of  Gar  Vc,a  complex  formations  by  other  approaches  such 
as  Monte  Carlo  method. 


Figure  9:  Calculated  energy  profile  of  the  gallium  Frenkel  defects.  Starting  from  the  perfect 
GaN  bulk  to  the  1st  (solid  dots)  and  2nd  Ga-tca  (solid  triangles)  configurations. 
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CN  and  SiN  complexes  in  dilute  nitrides  (GaAs:N  and  GaP:N) 

In  dilute  nitrides,  the  nitrogen  atoms  are  not  fitted  well  into  the  anion  site  due  to 
its  small  size.  As  a  result,  the  nitrogen  is  highly  strained.  On  the  other  hand,  the 
diatomic  molecule  such  as  N2,  NO,  O2,  CN  or  SiN  is  slightly  bigger  and  can  substitute  on 
the  anion  size  with  a  reasonable  energy.  As  a  result,  when  there  are  impurities  or 
intentionally  doped  elements,  such  as  C  or  Si,  the  diatomic  substitution  on  the  anion  site 
such  as  CNas  and  SiNAs  can  form. 

CN  complexes  in  diluted  nitrides 

CN  molecules  in  GaAs  and  GaPhave  been  studied.  Two  possible  lattice 
locations:  (1)  CN  molecules  substituting  for  anions  (CNas  or  CNp),  as  shown  in  Fig.  10 
and  (2)  interstitial  CN  molecules  at  the  Td  site  surrounded  by  four  Ga  atoms  (CN,)  are 
investigated.  All  possible  charge  states  and  various  orientations  were  considered. 
Calculated  formation  energies  show  that  the  molecule  favors  substituting  for  the  anion 
site  over  the  interstitial  configuration  for  all  equilibrium  growth  conditions  with  a  margin 
of  at  least  1  eV.  The  calculations  predict  the  CN  molecule  to  produce  a  level  with  a 
strong  ppn *  molecular  orbital  characteristic  at  approximately  0.5  -  0.8  eV  above  the 
VBM.  In p- type  conditions,  where  the  Fermi  level  is  located  below  this  ppn*  level,  the 
molecules  are  triply  bonded  and  form  donor  defect  centers  (double  donor  CN^  and 
single  donor  CNj.+ ).  The  calculations  with  full  relaxation  show  that  triply  bonded  CN^S 
and  CNj+  do  not  orient  symmetrically  in  the  Ga  tetrahedron  but  instead  are  tilted  in  order 
to  gain  better  interactions  with  neighboring  Ga  atoms.  The  calculated  vibration 
frequencies  of  the  triply  boned  CN^  and  CN;I+  are  in  reasonable  agreement  with  the 

measurement  by  Ulrici  and  Clerjaud  [21]  supporting  their  identification  that  the  bonding 
is  triple-bond  type.  Although  neither  CNas  nor  CN,  have  the  C-N  bond  oriented  along  the 
<100>  directions  as  proposed  by  Ulrici  and  Clerjaud,  the  low  rotation  barrier  of  CN^ 
makes  it  possible  that  the  molecule  might  be  constantly  rotating,  leading  to  an  average 
orientation  in  the  <100>  direction.  At  higher  Fermi  level,  the  CN)  can  accept  an 

electron  and  becomes  CN”  whereas  CN^  can  accept  1,  2,  3,  or  4  electrons  (depending 
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on  the  Fermi  level)  and  becomes  CN'a+s  ,  CN”S  (metastable),  CNAs  ,  or  CNAs , 

respectively.  We  found  that  the  C-N  bond  length  and  vibrational  frequency  change 
linearly  (to  a  good  approximation)  with  the  number  of  electrons  added  into  the 
antibonding  states  of  the  molecule. 

Further  details  of  this  work  can  be  found  in  an  attached  manuscript  [P2], 


O 

O 


Figure  10:  Local  atomic  geometry  of  a  CN  molecule  substituting  on  an  As  site  in  GaAs:  (a) 
symmetric  configuration,  CNAs(sym)  and  (b)  asymmetric  configuration,  CNAs(asym).  Both  are 
shown  for  charge  state  2+. 


SiN  complexes  in  diluted  nitrides 

Yu  et  al.  [22]  reported  that  when  GaAsi_xNx  is  doped  with  Si,  the  interaction 
between  Si  and  N  leads  to  a  “mutual  passivation,”  characterized  by  an  increase  in  the 
band  gap  and  an  elimination  of  the  electrical  activity  of  the  Si  donor.  We  found  that  the 
SiN  molecule  substituting  on  the  As  site  (SiNAs)  provide  a  better  model  to  explain  the 
mutual  passivation  of  shallow  donor  Si  and  isovalent  N  in  dilute  GaAsN  alloys.  Instead 
of  the  recently  proposed  pairing  of  Si  and  N  on  adjacent  substitutional  sites  (SiGa-NAs) 
[22,  23]  we  find  that  N  changes  the  behavior  of  Si  in  dilute  nitride  alloys  in  a  more 
dramatic  way.  N  and  Si  combine  into  a  deep-acceptor  split  interstitial,  where  Si  and  N 
share  an  As  site  SiNAs,  with  a  significantly  lower  formation  energy  than  that  of  the  Sioa- 
Nas  pair  in  n-type  GaAs  and  dilute  GaAsN  alloys.  The  formation  of  SiNAs  explains  the 
GaAs  band-gap  recovery  and  the  appearance  of  a  photoluminescence  peak  at  0.8  eV. 
This  model  can  also  be  extended  to  Ge-doped  GaAsN  alloys,  and  correctly  predicts  the 
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absence  of  mutual  passivation  in  the  case  of  column- VI  dopants.  This  work  has  been 
published  in  Phys.  Rev.  Lett.  [12] 


Figure  11:  Local  structure  of  (a)  (SiGa-NAs)+  and  (b)  (SiNAs)  .  In  the  SiGa-NAs  configuration, 
SiGa  is  a  nearest  neighbor  of  NAs;  in  the  SiNAs  configuration,  Si  and  N  share  an  As  site. 


Local  structure  of  Indium  oxynitride  (InON)  from  x-ray  absorption  spectroscopy 

Indium  oxynitride  can  be  considered  as  an  alloy  between  indium  nitride  (InN)  and 
indium  oxide  (h^Cb).  The  two  parent  compounds  have  very  different  crystal  structures. 
In  this  work,  the  local  structures  of  indium  oxynitride  (InOxNi.x)  with  varied  O  contents, 
prepared  by  rf  magnetron  sputtering,  were  identified  by  using  a  combination  of  first- 
principles  calculations  and  synchrotron  x-ray  absorption  near  edge  structures  (XANES). 
The  InOxNi.x  samples  were  characterized  by  indium  L3-edge  XANES  measurements  at 
Synchrotron  Light  Research  Institute  in  Thailand.  To  understand  the  local  structures  of 
the  samples,  we  independently  performed  first-principles  indium  L3-edge  simulations  of 
wurtzite  InN,  wurtzite  InNo.60o.4,  and  bixbyite  R12O3  crystals.  The  measured  XANES 
spectrum  of  sample,  with  no  intentionally  doped  O,  is  in  agreement  with  the  simulated 
spectrum  of  pure  wurtzite  InN.  The  spectrum  of  the  sample,  which  contains  only  O,  is 
also  in  agreement  with  the  simulated  spectrum  of  pure  bixbyite  h^Cb.  Interestingly,  the 
spectrum  of  sample  containing  ~40%  of  O  (InNo.60o.4)  is  almost  identical  to  that  of  pure 
InN  except  at  low  energy  range.  This  provides  strong  evidence  that  the  structure  of  the 
sample  remains  wurzite  like  that  of  pure  InN.  For  the  sample  containing  more  than  40% 
of  O,  the  spectra  show  a  mixed  signature  between  that  of  wurtzite  (InNo.60o.4)  and 
bixbyite  (I^Cf)  structures.  This  work  illustrates  that  XANES  is  a  powerful  technique  for 
identification  of  local  structures,  especially  when  it  is  used  in  conjunction  with  first 
principles  calculations. 

Further  details  of  this  work  can  be  found  in  an  attached  manuscript  [PI]. 
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Hydrogen  site  in  ZnO:  the  role  of  Ca  impurity 

In  2000,  Van  de  Walle,  proposed  that  H  prefers  to  form  a  strong  O-H  bond  with 
lattice  O  and  acts  exclusively  as  a  donor,  H+  [24],  This  behavior  is  different  from  H  in 
almost  all  other  semiconductors  in  which  H  is  amphoteric,  counteracting  the  dominant 
dopant  [25].  The  prediction  that  H  is  an  exclusive  donor  stirred  up  considerable  research 
activities,  because  it  may  explain  the  difficulty  in  the  /^-doping  [26-34],  Although,  it  is 
clear  that  H+  prefers  to  bind  strongly  with  O  at  a  bond  length  of  approximately  1  A  and 
with  a  vibration  frequency  around  3500  cm'1,  the  details  is  completely  controversial: 
McCluskey  et  al.  [29]  observed  an  IR  peak  at  3326  cm'1  in  samples  provided  by  Cermet, 
while  Lavrov  et  al.  [35]  observed  a  different  peak  at  3611  cm"1  in  samples  by  Eagle- 
Picher.  This  indicates  that  there  exist  two  forms  of  O-H  in  these  two  types  of  samples. 
These  3326  and  361 1 -cm'1  IR  peaks  are  consistent  with  the  calculated  frequencies  of  H  at 
the  antibonding  and  bond-center  sites,  respectively.  Based  on  first  principles 
calculations,  we  showed  that  Ca,  found  only  in  Cermet  samples,  plays  a  pivotal  role  in 
switching  the  H  from  the  bond-center  site  to  the  anti-bonding  site.  We  found  that  with 
Ca,  the  anti-bonding  site  is  0.5  eV  more  stable  than  the  ground-state  bond-center  site  in 
pure,  Ca-free  ZnO. 

Further  details  of  this  work  can  be  found  in  an  attached  manuscripts  [P3][P4], 
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Figure  13:  Schematic  drawings  the 
four  H+  sites  in  Ca-free  ZnO  and 
near  a  Ca  impurity.  Calculated 
atomic  relaxations  of  Ca-AB0±  (c)  in 
the  absence  and  (d)  in  the  presence 
of  H+. 
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Examples  of  Future  Research  Activities 

Binding  between  MgAi  and  Mg int  (as  well  as  between  BeAi  and  Beint)  in  AIN 
In  order  for  Mg  and  Be  to  be  a  good  p- type  dopant  in  AIN,  it  is  also  important  that  the 
self  compensations  can  be  controlled.  In  GaN,  our  calculations  have  shown  that  Be  could 
be  a  superior  p-type  dopant  over  Mg  when  Be  is  substituting  for  Ga.  However,  due  to  a 
smaller  size  of  Be,  the  interstitial  form  of  Be  in  GaN  also  have  rather  low  formation 
energy,  especially  when  the  sample  is  p- type.  The  similar  study  will  be  repeated  for  AIN. 

Roles  of  native  defects  in  III-N  alloys 

In  this  program  the  native  defects  in  AIN  has  been  studied.  This  new  result,  combining 
with  the  earlier  studied  of  native  defects  in  GaN,  gives  us  a  more  complete  picture  of  the 
native  defects  in  this  class  of  compounds.  However,  most  Ill-nitride  devices  use  the 
alloys  not  just  the  parent  compounds.  Therefore,  native  defects  in  alloys  should  be 
directly  studied.  We  will  directly  study  the  properties  of  native  defects  in  alloys  in  the 
next  program. 

Statistical  analysis  of  Frenkel  pair  stability 

Based  on  the  calculated  energy  of  (Ga)  Frenkel  pair  vs  separation  of  the  pair  obtained 
from  this  program,  the  stability  of  the  Frenkel  pair  can  be  further  modeled  by  using 
Monte  Carlo  method.  From  the  energy  profile,  the  realistic  simulation  of  the  creation  of 
the  Frenkel  pair  by  electron  irradiation  under  different  flux,  electron  energy,  and 
temperature  of  the  system  can  be  made. 

Nitrogen  Frenkel  pair  in  GaN 

To  complete  the  Frenkel  pair  study,  the  nitrogen  Frenkel  pair  will  be  investigated. 
Nitrogen  Frenkel  pair  composed  of  a  nitrogen  interstitial  and  a  nitrogen  vacancy  forming 
a  pair.  The  properties  and  the  stability  of  nitrogen  Frenkel  pair  could  be  very  different 
from  Ga  Frenkel.  This  is  because  the  nitrogen  interstitial  is  actually  in  the  split 
interstitial  form  (two  N  atoms  forming  an  N2  and  sharing  a  nitrogen  site). 
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Synchrotron  x-ray  absorption  near  edge  structures  (XANES)  measurements  of  In  L3  edge  is  used  in 
conjunction  with  first  principles  calculations  to  characterize  rf  magnetron  sputtered  indium 
oxynitride  at  different  O  contents.  Good  agreement  between  the  measured  and  the  independently 
calculated  spectra  are  obtained.  Calculations  show  that  the  XANES  spectra  of  this  ahoy  are  sensitive 
to  the  coordination  numbers  of  the  In  atoms,  i.e.,  fourfold  for  indium  nitride-like  structures  and 
sixfold  for  indium  oxide-like  structures,  but  not  to  the  substitution  of  nearest  neighbor  N  by  O  or 
vice  versa.  ©  2008  American  Institute  of  Physics.  [DOI:  10.1063/1.2965802] 


Indium  oxynitride  can  be  considered  as  an  ahoy  between 
indium  nitride  (InN)  and  indium  oxide  (ln203).  This  is  unlike 
most  of  the  traditional  semiconductor  alloys  in  the  sense  that 
the  crystal  structures  as  well  as  the  anion  valencies  of  the 
two  parent  compounds  are  different.  InN  is  one  of  the  highly 
studied  Ill-nitride  semiconductors  (GaN,  AIN,  and  InN)  be¬ 
cause  Ill-nitride  alloys  are  widely  used  for  optoelectronic 
applications.  Despite  real  applications  and  substantial  re¬ 
search  in  Ill-nitride  materials,  InN  bandgap  has  been  mis¬ 
taken  to  be  —1.9  eV  for  a  long  time  (see  the  discussion  in 
Refs.  1  and  2).  Only  recently,  the  actual  InN  bandgap  of 
—0.7  eV  has  been  realized.3,4  The  apparently  large  bandgap 
observed  in  the  past  is  most  likely  due  to  the  Moss-Burstein 
shift  caused  by  substantial  unintentional  carriers  that  are 
typical  for  InN  grown  by  traditional  techniques  (dc  discharge 
or  sputtering).1,5-7  Because  InN  has  a  (now-realized)  small 
bandgap  of  —0.7  eV  and  ln203  has  a  large  optical  gap  of 
—  3.6  eV,8  the  band  gap  of  indium  oxynitride  can  potentially 
be  engineered  in  a  very  wide  range.  This,  combining  with  the 
availability  of  low  temperature  growth  techniques,  such  as  rf 
magnetron  sputtering,9  make  indium  oxynitride  a  strong  can¬ 
didate  for  optical  coating  applications. 

For  alloys  with  low  O  content,  it  has  been  reported  that 
O  atoms  substitute  for  N  in  the  wurtzite  InN  crystal 
structure.10  However,  in  higher  O  content  alloys  the  crystal 
structure  remains  unclear.  Previously,  an  attempt  to  study 
indium  oxynitride  structures  by  N  X-edge  x-ray  absorption 
near  edge  structures  (XANES)  has  been  done.11  However, 
probing  local  structure  of  anions  in  this  alloy  system  is  not 
the  most  direct  way  because  anions  in  both  InN  and  ln203 
are  fourfold  coordinated.  Unlike  anion,  the  coordination 
number  of  In  in  the  two  compounds  are  different,  i.e.,  it  is 
sixfold  in  ln203  [Fig.  1(a)]  and  fourfold  in  InN  [Fig.  1(b)]. 
Therefore,  probing  the  local  structure  of  In  atoms  should 
give  more  direct  information  on  the  ahoy  structures.  In  this 
letter,  the  In  L3-edge  XANES  is  used  to  characterize  indium 
oxynitride  with  varied  compositions.  The  XANES  has  been 
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proven  to  be  a  powerful  tool  in  resolving  the  local  structures 
around  the  absorbing  atoms  (in  this  case  In  atoms).12-14 

Indium  oxynitride  films  were  grown  by  rf  magnetron 
sputtering  (Edwards  Auto  306)  at  room  temperature  using  a 
technique  called  reactive  gas  timing.15,16  The  growth  process 
started  with  the  preevacuation  of  the  chamber  to  the  order  of 
10-5  Pa.  Then  N2  and  02  gas  were  flown  interchangeably  at 
the  flow  rate  of  10  standard  cubic  centimeter  per  minute  onto 
the  99.999%  purity  In  target.  The  sputtering  gas  pressures 
were  set  at  0.34  and  0.32  Pa  for  N2  and  02,  respectively.  By 
controlling  N2  and  02  gas  timing,  indium  oxynitride  samples 
with  varied  O  contents  were  obtained.  The  interval  of  gas 
timing  cycles  used  for  each  sample  is  shown  in  Table  I  and 
the  samples  are  named  according  to  the  gas  timing  cycles. 
The  rf  plasma  power  was  set  at  100  W.  All  samples  were 
grown  on  polyethylene  terephthalate  substrate  to  the  thick¬ 
ness  of  —1  /zm  (typical  growth  time  is  approximately  1  h). 
To  measure  the  actual  O  content  in  each  sample,  the  [0]:[ln] 
ratios  were  determined  using  Auger  electron  spectroscopy 
(AES)  with  capability  of  sputter  depth  profiling.  The  O  con¬ 
tents  were  found  to  be  homogeneous  throughout  the  him 
thickness.  Note  that  AES  measured  O  in  ah  forms.  The  op¬ 
tical  bandgaps  measured  by  UV-visible  spectroscopy  were 
found  to  increase  with  02  gas  timing  ratio,  as  shown  in  Table 
I.  This  is  consistent  with  the  work  of  Yi  et  al.n  that  reported 
the  increase  in  optical  bandgap  of  indium  oxynitride  with  O 
content.  Although  the  optical  bandgap  of  N30/O0  sample  is 


FIG.  1.  (Color  online)  The  local  structure  around  In  atoms  used  in  the 
simulations  of  the  In  L3-edge  XANES  of  (a)  ln203  and  (b)  InN.  All  bond 
distances  are  given  as  a  percentage  difference  from  an  average  ln203  bond 
distance  (dcalc=2.170  A). 
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TABLE  I.  The  gas  timing  conditions,  optical  bandgaps,  and  [0]:[ln]  com¬ 
positions  of  the  samples. 

Gas  timing  (s) 

Optical 
bandgap  (eV) 

[O] :  [In] 
composition 

Sample 

n2 

o2 

N30/O0 

30 

0 

1.5 

0.16 

N30/O5 

30 

5 

1.6 

0.43 

N30/O10 

30 

10 

1.6 

0.96 

N30/020 

30 

20 

2.4 

1.08 

N0/O30 

0 

30 

3.4 

1.50 

larger  than  that  of  pure  InN  and  the  sample  is  expected  to 
contain  a  few  percent  of  unintentionally  doped  O,  this  gap 
widening  in  low  O  samples  is  most  likely  due  to  Moss- 
Burstein  effects1,5-7  not  the  changes  in  the  crystal  structure. 
As  will  be  discussed,  L3-edge  XANES  of  In  is  sensitive  to 
the  change  in  In  coordination  number  but  not  to  the  inter¬ 
changing  between  N  and  O.  Therefore,  XANES  of  InN  with 
a  few  percent  of  N  replaced  by  O  are  expected  to  be  almost 
identical  to  pure  InN  because  they  are  structurally  the  same. 

The  samples  were  characterized  by  In  L3-edge  XANES 
measurements  in  the  fluorescent  mode  with  a  13 -component 
Ge  detector  (Canbera)  at  the  x-ray  absorption  spectroscopy 
beamline  (BL-8)  of  the  Siam  Photon  Source  (electron  energy 
of  1.2  GeV,  beam  current  120-80  mA),  National  Synchro¬ 
tron  Research  Center,  Thailand.  Double  crystal  monochro¬ 
mator  Si  (111)  was  used  to  scan  the  synchrotron  x  ray  with 
the  photon  energy  step  of  0.25  eV  in  the  range  of 
3700  to  3850  eV,  covering  the  XANES  region  of  In  L3  edge. 
The  measured  spectra  are  shown  in  Fig.  2(a). 

In  order  to  understand  the  local  microscopic  structure, 
we  performed  first  principles  In  L3-edge  XANES  simulations 
of  wurtzite  InN  and  bixbyite  ln203  crystals.  The  detailed 
crystal  structures  were  optimized  based  on  first  principles 
pseudopotential  calculations.  We  used  density  functional 
theory  with  local  density  approximation  and  ultrasoft 
pseudopotentials  as  implemented  in  the  VASP  code.17  The 
cutoff  energy  for  the  plane  wave  basis  set  was  set  at  400  eV. 


FIG.  2.  (Color  online)  (a)  Normalized  In  L3-edge  XANES  spectra  of  indium 
oxynitride  samples  prepared  under  different  gas  timing  with  increased  02 
timing  from  bottom  curve  to  top  curve.  Circles  show  recorded  data  points 
and  the  curves  are  the  (noise  removed)  fit  to  the  data.  Inset:  the  plots  from 
the  main  panel  without  offset  at  the  energy  range  near  the  shoulders  to 
highlight  the  changes  in  the  shoulder  height,  (b)  The  calculated  In  L3-edge 
XANES  spectra  of  wurtzite  InN,  bixbyite  ln203,  wurtzite  InN06O04,  and  the 
simulated  separated  phase  alloys  with  38%  and  49%  fourfold  In  atoms  (see 
text  for  detail). 


We  used  Monkhorst-Pack  k-point  mesh  for  Brillouin  zone 
integration  (7X7X7  for  InN  and  3X3X3  for  ln203).  All 
atoms  were  allowed  to  relax  until  the  residue  forces  were  less 
than  10-3  eV/A.  The  fully  relaxed  local  structures  surround¬ 
ing  In  atoms  are  shown  in  Figs.  1(a)  and  1(b)  for  ln203  and 
InN,  respectively.  All  In  atoms  in  wurtzite  InN  crystal  are 
equivalent.  However,  there  are  two  species  of  In  atoms  in 
ln203,  labeled  Ini  and  In2,  with  the  composition  ratio  of  1:3. 
These  structural  relaxations  by  first  principles  calculations  is 
a  crucial  step  to  obtain  the  computationally  unstrained  struc¬ 
tures  to  be  further  used  in  XANES  calculations. 

To  calculate  ab  initio  XANES  of  InN  and  ln203  based 
on  the  relaxed  crystal  structures  above,  we  used  FEFF8.2 
codes.18,19  The  codes  utilize  a  full  multiple  scattering  ap¬ 
proach  based  on  ab  initio  overlapping  muffin-tin  potentials. 
The  muffin-tin  potentials  were  obtained  using  self-consistent 
calculations  with  Hedin-Lundqvist  exchange-correlation 
function.  The  self-consistent  calculations  were  performed  in 
the  sphere  radius  4  A  (—40  atoms)  around  the  absorber  In 
atom,  which  is  more  than  sufficient  to  allow  charge  to  fluc¬ 
tuate  yielding  realistic  electron  relaxations.  The  full  multiple 
scattering  calculations  include  all  possible  paths  within  a 
larger  cluster  radius  of  7.4  A  (—140  atoms). 

The  calculated  In  L3-edge  spectra  of  InN  and  ln203  crys¬ 
tals  are  shown  in  Fig.  2(b)  (bottom  and  top  curves,  respec¬ 
tively).  For  ln203,  because  there  are  two  species  of  non¬ 
equivalent  In  atoms,  the  spectrum  shown  is  the  weight 
averaged  between  the  spectra  from  the  two  species.  Note  that 
the  spectra  obtained  from  the  two  species  of  In  in  ln203  are 
very  similar.  The  spectra  of  fourfold  In  in  InN  and  sixfold  In 
in  ln203  are,  however,  clearly  different;  making  In  L3-edge 
XANES  a  suitable  tool  to  study  the  local  structure  around  In 
atoms.  Our  test  calculations  also  show  that  the  XANES  spec¬ 
tra  of  the  fourfold  (sixfold)  In  atom  is  not  very  sensitive  to 
the  substitution  of  a  neighboring  N  by  O  (O  by  N).20 

The  calculated  spectra  of  InN  and  ln203  are  in  full 
agreement  with  the  measured  spectra  from  N30/O0  and  NO/ 
030  samples,  respectively.  The  main  features  of  InN 
XANES  spectrum  is  composed  of  two  shoulders  and  three 
peaks  labeled  as  s2,  pi,  p2,  and  p3.  All  five  features  are  in 
a  very  good  agreement  with  the  corresponding  features  in  the 
N30/O0  sample  in  both  positions  and  shapes.  Similar  agree¬ 
ment  can  be  found  between  the  spectrum  of  N0/O30  sample 
and  the  calculated  ln203  spectrum.  This  suggests  that  struc¬ 
turally  In  atoms  in  N30/O0  are  mostly  fourfold  (as  those  in 
InN)  and  in  N0/O30  are  mostly  sixfold  (as  those  in  ln203).  If 
we  compare  the  measured  spectra  of  different  O  contents 
[from  the  bottom  curve  to  the  top  curve  of  Fig.  2(a)],  we  can 
see  that  the  five  features  are  progressively  evolved  as  O  con¬ 
tent  increased.  The  progressive  changes  include:  (1)  the  re¬ 
duction  in  Si  peak  (for  clarity,  the  same  set  of  spectra  without 
offset  is  shown  in  the  inset),  (2)  the  broadening  and  increase 
in  magnitude  of  px  peak,  (3)  the  reduction  in  p2  peak,  and  (4) 
the  significant  broadening  accompanied  by  a  slight  shift  to 
higher  energy  of  peak  p3.  Interestingly,  the  spectrum  of 
N30/O5  (at  43%  O  content)  is  still  almost  perfectly  over¬ 
lapped  with  that  of  N30/O0  with  the  exception  of  sx  shoul¬ 
der,  which  is  slightly  reduced  [see  inset  of  Fig.  2(a)].  This 
indicates  that  most  of  In  atoms  in  N30/O5  remain  fourfold 
and  O  atoms  are  substituting  on  the  N  sites.  To  test  this 
assumption,  we  calculated  the  XANES  spectra  for  InN0  6O0  4 
alloy  in  the  wurtzite  structure  [the  second  curve  from  bottom 
in  Fig.  2(b)].  The  calculated  curve  is  almost  overlapped  with 
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that  of  (calculated)  pure  InN  for  the  photon  energy  beyond 
the  Si  shoulder.  Moreover,  a  small  reduction  in  sx  shoulder  of 
the  N30/O5  spectrum  in  comparison  with  that  of  N30/O0  is 
nicely  reproduced  by  the  InN06O04  alloy  simulation.  The 
reduction  in  the  sx  shoulder  can  be  attributed  to  the  changes 
in  the  lowest  conduction  band  states  due  to  the  replacement 
of  N  by  O  atoms.  Note  that  a  calculated  spectrum  based  on  a 
combination  of  fourfold  and  sixfold  In  (not  shown)  does  not 
give  such  a  good  agreement.  For  samples  with  higher  O  con¬ 
tents  (N30/O10  and  N30/020)  the  spectra  show  a  mixed  sig¬ 
nature  between  that  of  fourfold  and  sixfold  In  atoms.  The 
excess  O  atoms  above  the  substitutional  solubility  limit  (as¬ 
sumed  to  be  —40%)  of  InN  can  either  form  a  phase  sepa¬ 
rated  ln203-like  structure  or  an  inclusion  of  local  ln203-like 
structure  both  of  which  leads  to  a  formation  of  sixfold  In 
atoms.  In  this  model,  the  composition  in  the  sample  can  be 
written  as  (lnN0  6O0A)F(lnOl  5)l_F,  where  F  is  a  fraction  of 
fourfold  In  atoms  out  of  all  In  atoms.  The  value  of  F  can  be 
calculated  using  a  relationship  [O] :  [In]  =  0AF+ 1.5(1  -F), 
where  the  [0]:[ln]  ratio  of  each  sample  was  determined  from 
AES  and  is  shown  in  Table  I.  This  gives  the  F  values  of  0.49 
and  0.38  for  N30/020  and  N30/O10,  respectively.  The  cor¬ 
responding  spectra,  which  are  calculated  using  weight  aver¬ 
aged  between  InN06O04  and  pure  ln203,  are  shown  in  Fig. 
2(b).  The  spectra  show  the  mixed  signature  of  InN  and  ln203 
that  are  consistent  with  the  measured  spectra  of  N30/O10 
and  N30/020.  sx  shoulder  also  further  decreased  from  pure 
InN  as  more  O  content  is  added  (insets  of  Fig.  2). 

In  a  detail  investigation,  the  spectra  of  both  N30/O10 
and  N30/020  samples  appear  to  be  very  similar  to  each  other 
and  the  sx  shoulder  of  N30/020  is  only  slightly  lower  than 
N30/O10  [see  the  inset  of  Fig.  2(a)].  The  simulations  also 
show  the  same  trend,  i.e.,  the  sx  shoulder  of  the  spectrum 
with  F=0.38  is  just  slightly  lower  than  that  of  F=0.49  [see 
the  inset  of  Fig.  2(b)].  However,  when  observing  the  de¬ 
crease  in  Si  shoulder  in  comparison  to  that  of  InN,  we  can 
see  that  the  simulations  (F=0.38  and  F=0.49  spectra)  give 
larger  drop  in  the  sx  shoulder  in  comparison  to  those  of 
samples  N30/O10  and  N30/020.  This  suggests  that  N30/O10 
and  N30/020  samples  may  contain  more  fourfold  In  atoms, 
i.e.,  higher  F  ratios,  than  those  used  in  the  simulations.  We 
found  that  a  simulation  with  F  — 2/3,  i.e.,  2/3  of  In  atoms 
are  fourfold  and  1/3  of  In  atoms  are  sixfold,  can  provide  a 
better  agreement  in  the  S\  feature  with  the  experiment  (not 
shown).  We  speculate  that  there  exist  an  intermediate  crystal 
structure  that  have  a  fix  number  of  fourfold  In  atoms  and 
sixfold  In  atoms,  in  the  unit  cell.  Although  N  atoms  would 
prefer  to  decorate  fourfold  In  atoms  whereas  O  atoms  would 
prefer  sixfold  In  atoms,  the  interchanging  of  the  anions  is 
possible,  leading  to  a  very  different  optical  properties  in  the 
two  samples.  Since  the  XANES  signatures  are  not  very  sen¬ 
sitive  to  these  anion  interchanging,  the  spectra  of  the  two 
samples  remain  very  close.  Further  work  is  needed  to  settle 
down  this  issue. 

In  summary,  we  have  performed  In  L3-edge  XANES 
study  of  indium  oxynitride  prepared  by  rf  magnetron  sputter¬ 
ing  with  different  oxygen  contents.  The  independent  first 
principles  XANES  calculations  of  bulk  InN  and  ln203  give 
unambiguous  agreement  with  the  spectra  from  the  samples 
prepared  by  using  only  N2  and  02  gas,  respectively.  Our 
results  suggest  that  as  much  as  40%  of  oxygen  can  replace  N 
in  fourfold  InN  structure.  The  spectra  from  the  samples  with 


higher  O  contents  show  that  the  alloy  contains  both  fourfold 
and  sixfold  indium  atoms.  This  implies  that  either  the 
samples  are  phase  separated  or  there  exist  a  preferred  inter¬ 
mediate  crystal  structure  containing  both  fourfold  and  sixfold 
In  atoms.  Our  calculations  also  show  that  the  XANES  fea¬ 
tures  are  not  strongly  affected  by  the  substitution  of  N  by  O 
or  vice  versa. 
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A  carbon-nitrogen  molecule  in  GaAs  and  GaP  is  investigated  by  using  first-principles  density  functional 
pseudopotential  calculations.  The  formation  energy  calculations  show  that  the  molecule  favors  substituting  for 
an  anion  site  (As  or  P)  over  being  an  interstitial  under  all  equilibrium  growth  conditions.  Under  p-type 
conditions,  the  molecule  exhibits  the  characteristic  triple  bonding  and  acts  as  a  donor.  When  the  Fermi  level  is 
higher  ( n  type),  the  molecule  acts  as  an  acceptor  by  accepting  electrons  into  its  antibonding  states  and  exhibits 
double  or  single  bonding.  The  bond  length  and  vibrational  frequencies  for  each  configuration  are  calculated 
and  compared  to  those  in  recent  experiments.  Trends  in  the  changes  in  bond  length  and  vibrational  frequencies 
with  respect  to  the  number  of  electrons  in  the  antibonding  states  are  discussed. 
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I.  INTRODUCTION 

Small  amounts  of  nitrogen  can  strongly  affect  the  elec¬ 
tronic  properties,  such  as  the  band  gap  and  electron  effective 
mass,  of  GaP  and  GaAs.  This  opens  up  an  opportunity  to 
continuously  modify  the  band  gap  (as  well  as  other  proper¬ 
ties)  of  these  III-V  semiconductors  for  advanced  electronic 
and  optoelectronic  devices.  This  prospect  has  stimulated  sub¬ 
stantial  research  on  dilute  GaAsN  and  GaPN  alloys  both 
experimentally1-7  and  theoretically.8-14  It  is  therefore  very 
important  that  defects  and  impurities  in  these  alloys  are  stud¬ 
ied  in  detail.  Carbon  is  a  common  impurity  in  various  growth 
techniques15,16  and  sometimes,  it  is  intentionally  used  as  a 
/7-type  dopant;17-20  complex  formations  between  C  and  N 
could  therefore  be  an  important  issue.  Recently,  Ulrici  and 
Clerjaud21  observed  a  sharp  local  vibration  mode  at 
2087  cm-1  (at  T=T  K)  in  GaP  and  a  similar  mode  in  GaAs. 
They  identified  it  as  a  CN  complex  with  a  triple  bond  (be¬ 
tween  C  and  N)  aligned  along  the  [100]  direction. 

An  experimental  work21  rigorously  identified  the  chemi¬ 
cal  composition  and  the  dipole  direction  of  the  complex; 
however,  the  actual  location  of  the  CN  molecule  in  the  lattice 
is  still  unclear.  Both  CN  substituting  on  anion  sites  (CNAs  in 
GaAs  or  CNP  in  GaP)  and  interstitial  CN  (CN;)  are  poten¬ 
tially  consistent  with  experiment  and  were  proposed  as  pos¬ 
sible  models  in  Ref.  21.  In  other  semiconductor  compounds 
(such  as  ZnO,  GaN,  and  ZnSe),  small  diatomic  molecules 
such  as  N2,  02,  or  CN  have  been  studied  both  by  first- 
principles  calculations  and  by  different  experiments  and 
found  to  prefer  substitution  on  an  anion  site.15,22-25  However, 
the  energy  difference  between  the  substitutional  and  intersti¬ 
tial  sites  varies  from  material  to  material;  explicit  calcula¬ 
tions  are  therefore  required  to  determine  which  site  is  more 
favorable  for  a  given  material. 

In  this  paper,  we  use  first-principles  calculations  to  calcu¬ 
late  the  formation  energies  of  CNAs(orP)  and  CN;.  We  find 
that,  similar  to  other  semiconductor  compounds  that  have 
been  studied  before,  CN  energetically  prefers  substituting  on 
anion  sites  of  GaAs  and  GaP  over  interstitial  sites  in  all 
thermal  equilibrium  growth  conditions.  We  also  calculate  the 
vibrational  frequency  of  the  complex  and  find  reasonable 


agreement  with  the  measured  value.  The  CN  molecule  be¬ 
haves  very  similarly  in  GaAs  and  GaP,  and  therefore,  we  will 
focus  our  discussion  on  CN  in  GaAs.  Numerical  results  will 
be  reported  for  both  GaAs  and  GaP. 

II.  COMPUTATIONAL  METHOD 

Our  calculations  are  performed  by  using  density  func¬ 
tional  theory  with  the  local  density  approximation  (LDA).26 
We  use  ultrasoft  pseudopotentials,27  as  implemented  in  the 
VASP  code.28  The  Ga  3d  electrons  are  treated  as  valence  elec¬ 
trons.  The  cutoff  energy  for  the  plane- wave  basis  set  is  262 
eV.  We  use  a  supercell  with  64  atoms  for  the  defect  studies, 
and  a  2  X  2  X  2  shifted  Monkhorst-Pack  special  k-point  grid 
(r  point  not  included)  is  employed  for  the  Brillouin  zone 
integration.  The  calculated  GaAs  lattice  constant  (5.60  A, 
which  is  within  1%  of  the  experimental  value)  is  used.  The 
electronic  deep  levels  introduced  by  CN  are  examined  by 
taking  averages  over  the  special  k  points.  All  the  atoms  are 
relaxed  by  minimization  of  the  Hellmann-Feynman  forces 
until  all  the  forces  are  less  than  0.05  eV/A. 

The  formation  energies  of  CN;  and  CNAs  are  defined,  fol¬ 
lowing  Refs.  22  and  29,  as 

E/C N?)  =  £tot(CN?)  -  £tot(bulk)  -  Me  -  /%  +  qEF, 

£/(cnas)  =  £tot(CN^s)  -  £tot(bulk)  -  Me  -  Mn  +  Mas  +  qEP, 

(1) 

where  Etot(Dq)  is  the  total  energy  of  the  supercell  with  defect 
D  in  charge  state  q.  £tot(bulk)  is  the  total  energy  of  the  su¬ 
percell  without  a  defect.  /zc,  /zN,  and  julAs  are  the  chemical 
potentials  of  a  C  atom,  an  N  atom,  and  an  As  atom,  respec¬ 
tively.  For  convenience,  julc  and  /zAs  are  referenced  to  the 
energy  of  a  C  atom  in  diamond  and  an  As  atom  in  solid  As, 
respectively.  juln  is  referenced  to  the  energy  of  an  N  atom  in 
a  free  N2  molecule.  Growth  conditions  for  GaAs  are  close  to 
equilibrium,  requiring  that  /^Ga+yaAs-yaGaAs--0.69  eV  (cal¬ 
culated  GaAs  heat  of  formation),  where  /xGa  is  a  chemical 
potential  of  a  Ga  atom,  which  is  referenced  to  the  energy  of 
a  Ga  atom  in  bulk  Ga.  EF  is  the  electron  Fermi  level  with 
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respect  to  the  valence-band  maximum  (VBM). 

The  local  vibrational  mode  frequencies  of  the  molecules 
are  calculated  by  using  the  so-called  frozen  phonon 
approach,30  in  which  we  followed  a  practical  methodology 
described  in  Ref.  31.  Because  C  and  N  atoms  have  similar 
masses,  we  equally  displaced  C  and  N  to  calculate  the  po¬ 
tential  energy  of  stretching  and  compressing  bonds.  The 
mass  of  the  oscillator  is  taken  as  the  reduced  mass: 


A*'  ?  V-W 

Jfl] +  TYIq 

where  mc  and  raN  are  the  masses  of  C  and  N  atoms.  We  find 
that  the  (harmonic)  stretching  vibrational  frequencies  of  the 
oscillator  agree  very  well  with  a  full  dynamic  matrix  calcu¬ 
lation  with  the  matrix  constructed  by  displacing  each  and 
every  atom  in  the  supercell  one  at  a  time  in  all  three  principal 
axes  directions.  A  test  calculation  (on  CNA+)  shows  that  the 
agreement  between  the  full  dynamic  and  reduced-mass  re¬ 
sults  is  better  than  5  cm-1.  An  important  benefit  of  the 
reduced-mass  calculation  (besides  a  significant  reduction  in 
computational  effort)  is  that  it  allows  us  to  calculate  the  an- 
harmonic  part  of  the  vibration,  which  we  found  to  be  on  the 
order  of  20  cm-1.  The  potential  energy  is  evaluated  at  seven 
values  of  displacements,  with  a  maximum  amplitude  of  up  to 
10%  of  the  C-N  equilibrium  bond  distance.  To  test  the  accu¬ 
racy  of  the  approach,  we  calculated  the  stretch  frequency  of 
a  C-N  mode  in  the  HCN  molecule  and  obtained  the  value  of 
2057  cm-1,  which  is  slightly  smaller  than  the  experimental 
value  of  2089  cm-1.32  The  calculated  equilibrium  C-N  bond 
distance  is  fortuitously  equal  to  the  experimental  value  of 
1.156  A.  The  small  underestimation  in  vibrational  fre¬ 
quency,  which  is  by  32  cm-1,  is  expected  to  be  systematic. 
Therefore,  we  have  added  the  value  of  coER=32  cm-1  to  all 
of  our  calculated  frequencies  as  a  systematic  correction. 

III.  RESULTS  AND  DISCUSSION 

We  investigated  two  forms  of  CN  molecules  in  GaAs:  at 
the  substitutional  As  site  and  at  the  interstitial  site.  If  initially 
the  CN  molecule  is  symmetrically  placed  at  the  As  site  (in 
the  form  of  a  split-interstitial  configuration)  with  its  principal 
axis  aligned  along  the  [100]  direction,  due  to  the  symmetry, 
the  orientation  of  the  molecule  will  remain  unchanged  even 
when  relaxation  is  allowed.  The  resulting  symmetric  con¬ 
figuration,  CNAs(sym),  is  shown  in  Fig.  1(a).  However,  this 
orientation  is  not  the  lowest  energy  one  for  all  charge  states. 
Breaking  the  symmetry  causes  the  molecule  (in  some  charge 
states)  to  spontaneously  rotate,  without  any  barrier,  into  an 
asymmetric  configuration  [Fig.  1(b),  CNAs(asym)  in  the  2+ 
charge  state]. 

For  an  interstitial  CN,  there  are  two  possible  sites  for  the 
zinc-blende  crystal:  the  tetrahedral  interstitial  sites  sur¬ 
rounded  by  either  Ga  [7^(Ga)]  or  As  atoms  [7^(As)].  The 
calculations  show  that  CN*  energetically  prefers  the  7^(Ga) 
site  over  the  Td( As)  site.  The  T^Ga)  site  is  surrounded  by 
four  Ga  atoms  in  a  tetrahedral  configuration  similar  to  the 
local  structure  surrounding  a  substitutional  As  lattice  site. 
Again,  if  the  molecule  is  initially  placed  in  a  high-symmetry 
orientation,  CN;(sym),  the  orientation  remains  fixed  by  sym¬ 


FIG.  1.  (Color  online)  Local  atomic  geometry  of  a  CN  molecule 
substituting  on  an  As  site  in  GaAs:  (a)  symmetric  configuration, 
CNAs(sym),  and  (b)  asymmetric  configuration,  CNAs(asym).  Both 
are  shown  for  the  2+  charge  state. 

metry,  but  breaking  the  symmetry  can  lead  to  rotation  of  the 
molecule  into  a  lower-energy  asymmetric  configuration, 
CN;  (asym). 

A.  Substitutional  CN  molecules  (CNAs) 

The  substitutional  molecule  CNAs  introduces  a  deep  donor 
level  (at  approximately  0.25  eV)  in  the  GaAs  band  gap,  as 
shown  in  Fig.  2.  The  charge  density  (Fig.  3)  of  these  levels 
resembles  the  ppif  molecular  orbital  of  the  CN  molecule.22 
The  molecular  orbital  theory  tells  us  that  the  CN  molecule 
has  the  following  molecular  orbitals,  which  are  in  order  of 
increasing  energy:  ssa ,  sscr *,  ppir  (doublet),  ppc t,  ppjf 
(doublet),  and  ppa* .  In  the  neutral  charge  state,  the  CN  mol¬ 
ecule  on  an  As  site  has  12  valence  electrons:  9  from  the  CN 
itself  and  3  more  contributed  from  the  surrounding  Ga  at¬ 
oms.  These  12  electrons  occupy  the  CN  molecular  orbitals 
up  to  the  lower  ppir*  level  and  leave  the  higher  ppir*  level 
empty.  Because  the  lower  ppir*  level  is  located  deep  inside 
the  GaAs  band  gap,  CNAs  can  donate  up  to  two  electrons  and 
become  CNA+  when  the  Fermi  level  of  the  system  falls  below 
the  0.25  eV  donor  level.  The  bonding  between  the  C  and  N 
atoms  in  CNA+  is  a  triple  bond  (C  =  N),  with  a  bond  distance 
of  1.17-1.20  A  (depending  on  the  orientation),  which  is  in 
good  agreement  with  the  known  triple-bond  length  between 
C  and  N  in  HCN  of  1.156  A.33 


(a)  CNAs  (b)  CN. 


FIG.  2.  (Color  online)  Schematic  illustrations  of  the  ppir* 
single-particle  energy  levels  of  the  CN  molecule  in  the  GaAs  band 
gap.  (a)  Substituting  on  an  As  site  (CNAs)  and  (b)  CN  on  an  inter¬ 
stitial  site  (CN;).  The  solid  dots  show  the  electron  occupation  for 
the  neutral  charge  state. 
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FIG.  3.  (Color  online)  (a)  Charge  density  of  the  pp  tt*  states  of 
CNAs  that  are  located  in  the  GaAs  band  gap.  (b)  Contour  plot  of  one 
of  the  states  in  (a)  in  a  (001)  plane  cutting  through  the  CN  mol¬ 
ecule.  (c)  Charge  density  of  one  of  the  ppir  states  of  CNAs  that 
resides  below  the  valence-band  maximum,  (d)  The  corresponding 
contour  plot  of  (c). 

In  the  2+  charge  state,  the  ppir*  level  is  unoccupied; 
hence,  the  CNA+  has  little  interaction  with  the  surrounding 
neighbor.  The  molecule  can  therefore  rotate  quite  freely,  and 
it  indeed  moves  away  from  the  symmetric  configuration  to  a 
new  orientation  that  has  a  lower  energy.  Our  calculated  en¬ 
ergy  profile  for  the  molecule  rotating  from  the  [001]  to  the 
(near)  [110]  direction,  which  is  based  on  the  nudged  elastic 
band  (NEB)  method,34  is  shown  in  Fig.  4.  We  can  see  that 
starting  from  the  symmetric  [001]  orientation,  the  molecule 
can  rotate  without  any  barrier  to  a  new  orientation  in  which 
the  molecule’s  principal  axis  is  nearly  parallel  with  the  [110] 
direction,  lowering  its  energy  by  0.27  eV  in  the  process.  By 
investigating  the  charge  density  (not  shown),  we  found  that 
the  asymmetric  configuration  allows  better  interactions 


FIG.  4.  (Color  online)  Calculated  energy  of  CNA+  as  the  CN 
molecule  rotates  from  the  asymmetric  to  symmetric  and  back  to  an 
equivalent  asymmetric  configuration.  The  symmetric  configuration 
is  the  saddle  point  along  this  path. 


TABLE  I.  Calculated  bond  lengths  dc_ N,  stretching  frequencies 
co,  and  formation  energies  Ef  (for  the  Fermi  level  at  the  VBM)  of 
CN  molecules  in  GaAs  in  the  symmetric  (asymmetric  in  parenthe¬ 
ses)  configurations.  For  each  charge  state,  the  lower-energy  con¬ 
figuration  is  shown  in  boldface. 


Charge  state 
(?) 

^C-N 

(A) 

CO 

(cm-1) 

Ef  (eV) 

Ga  rich 

As  rich 

CN  on  the  As  site  (CNAs) 

2+ 

1.20 

1779 

2.42 

3.11 

(1.17) 

(2052) 

(2.15) 

(2.84) 

1+ 

1.25 

1481 

2.42 

3.11 

(1.22) 

(1734) 

(2.46) 

(3.15) 

0 

1.30 

1328 

2.58 

3.27 

(1.25) 

(1345) 

(3.08) 

(3.77) 

1- 

1.35 

1173 

3.13 

3.82 

(1.31) 

(1182) 

(3.34) 

(4.03) 

2- 

1.38 

1095 

4.03 

4.72 

(1.40) 

(1097) 

(3.82) 

(4.51) 

CN  on  the  interstitial  site  (CN;) 

1+ 

1.20 

1897 

4.41 

4.41 

(1.19) 

(1934) 

(3.85) 

(3.85) 

0 

(1.25) 

(1519) 

(4.37) 

(4.37) 

between  the  bonding  pp  states  of  CN  and  two  of  the  sur¬ 
rounding  Ga  atoms,  resulting  in  a  lower  formation  energy. 

Starting  from  CNA+,  when  the  Fermi  level  is  raised,  the 
defect  levels  characteristic  of  ppir*  molecular  orbitals  be¬ 
come  occupied  with  electrons,  and  the  charge  state  of  CNAs 
increases  from  2+  to  1+,  0,  1-,  and  2-.  Because  these  defect 
states  correspond  to  antibonding  orbitals  of  the  free  mol¬ 
ecule,  occupying  them  leads  to  an  increase  in  the  CN  bond 
length  (of  approximately  0.05  A  per  electron  added),  as 
shown  in  Table  I.  In  addition,  the  charge  density  plots  [Fig. 
3(a)  and  3(b)]  reveal  that  the  ppir*  states  hybridize  with  the 
dangling-bond  states  of  the  surrounding  Ga  atoms  in  the  di¬ 
rections  in  which  the  density  lobes  are  pointing.  For  instance, 
the  lower  two  lobes  in  Fig.  3(b)  clearly  show  the  bond  for¬ 
mation  between  the  C  atom  and  two  Ga  atoms,  whereas  the 
upper  two  lobes  still  maintain  the  molecular  ppir*  features. 
Because  of  the  possibility  of  bond  formation  with  the  neigh¬ 
boring  Ga  atoms, CNA+,  CNAs,  and  CNA“  all  favor  the  sym¬ 
metric  configuration.35 

CNa“,  on  the  other  hand,  favors  an  asymmetric  configu¬ 
ration  (with  a  different  orientation  from  that  of  CNA+).  Be¬ 
cause  the  ppTr*  states  are  fully  occupied,  the  molecule  tries 
to  optimize  the  Coulombic  interaction  between  these  ppir* 
states  and  its  four  Ga  neighbors.  Since  the  charge  distribution 
for  the  CN  molecule  is  not  symmetric  (the  N  atom  is  more 
electronegative  than  the  C  atom),  the  asymmetric  configura¬ 
tion  allows  the  N  side  of  the  molecule  to  bind  with  three  Ga 
atoms  instead  of  two,  making  the  configuration  more  stable. 

The  formation  energies  of  CNAs  are  shown  in  Table  I 
(Table  II  for  CN  in  GaP)  and  plots  as  a  function  of  Fermi 
level  are  shown  in  Fig.  5  for  two  growth  conditions:  Ga  rich 
(/zGa=0)  and  As  rich  (yu,As  =  0).  As  can  be  seen  from  Eq.  (1), 
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TABLE  II.  Calculated  bond  lengths  dc_ N,  stretching  frequencies 
co,  and  formation  energies  Ef  (for  the  Fermi  level  at  the  VBM)  of 
CN  molecules  in  GaP  in  the  symmetric  (asymmetric  in  parentheses) 
configurations.  For  each  charge  state,  the  lower-energy  configura¬ 
tion  is  shown  in  boldface. 


Charge  state 
(q) 

^C-N 

(A) 

CO 

(cm-1) 

Ef  (eV) 

Ga  rich 

Prich 

CN  on  the  P  site  (CNp) 

2+ 

1.20 

1734 

2.25 

3.13 

(1-18) 

(2028) 

(2.12) 

(3.00) 

1  + 

1.25 

1517 

2.40 

3.28 

0 

1.31 

1352 

2.75 

3.63 

1- 

1.36 

1197 

3.58 

4.46 

(1.33) 

(1283) 

(3.74) 

(4.62) 

2- 

1.39 

1113 

4.77 

5.65 

(1.41) 

(1105) 

(4.50) 

(5.38) 

CN  on  the  interstitial  site  (CN;) 

1  + 

1.21 

1875 

4.43 

4.43 

(1.19) 

(1933) 

(3.96) 

(3.96) 

0 

(1.26) 

(1511) 

(4.82) 

(4.82) 

the  slope  of  the  plot  indicates  the  charge  state  of  the  CN 
molecule.  In  our  plots,  only  the  lowest  energy  charge  state 
(at  a  given  Fermi  energy)  is  shown.  For  example, 
CNAs(asym)  is  stable  in  the  2+  charge  state  when  the  EF  of 
the  system  is  located  between  0.0  and  0.25  eV;  therefore,  the 
line  has  a  slope  of  2.  For  0.25  eV<£F<0.70  eV,  the  sym¬ 
metric  configuration  becomes  more  stable  and  the  neutral 
and  1-  (for  a  narrow  range  near  £^=0.70  eV)  charge  states 
are  the  lowest  in  energy.  At  higher  EF ,  the  asymmetric  con¬ 
figuration  with  the  2-  charge  state,  i.e.,  CNA“(asym),  be¬ 
comes  the  most  stable.36  Therefore,  the  stable  configurations 
in  order  of  increasing  EF  are  as  follows:  CNA+(asym), 
CNAs(sym),  CNA“(sym),  and  CNA“(asym).  CNA+  is  never 
stable.  Because  CN  is  substituting  on  the  As  site  to  form 
CNAs,  its  formation  energy  is  lower  in  the  Ga-rich  condition 
than  that  in  the  As-rich  condition.  The  plots  also  show  that 
the  formation  energy  of  CNAs  is  always  lower  than  that  of 
CN;  for  both  the  Ga-  and  As-rich  conditions  throughout  the 
entire  Fermi-energy  range,  with  an  energy  difference  of  at 
least  1  eV.  This  implies  that  CNAs  is  far  easier  to  form  than 
CN;. 

The  trends  in  bonding  discussed  above  are  reflected  in  our 
calculated  values  of  vibrational  frequencies,  which  show  that 
the  frequency  is  highest  for  CNA+  (2052  cm-1  for  an  asym¬ 
metric  configuration  and  1779  cm-1  for  a  hypothetical  sym- 


FIG.  5.  (Color  online)  Formation  energies  of  CN  molecules  in 
GaAs  as  a  function  of  the  Fermi  energy.  The  solid  lines  and  dashed 
lines  represent  the  asymmetric  and  symmetric  configurations, 
respectively. 

metric  configuration).  This  is  because  in  the  2+  charge  state, 
the  bonding  states  of  the  molecule  are  all  occupied  and  the 
antibonding  states  ( ppir *)  are  all  empty.  As  electrons  are 
added  into  the  antibonding  states  of  CNAs,  the  C-N  bond  is 
weakened  and  the  frequency  is  reduced  to  1481,  1328,  1173, 
and  1097  cm-1  for  CNA+,  CNAs,  CNA“,  and  CNA“,  respec¬ 
tively. 

B.  Interstitial  CN  molecules  (CN;) 

The  CN  molecule  has  a  total  of  nine  electrons  filling  its 
molecular  orbitals  in  the  following  order:  ssa ,  ssa *,  ppir 
(doublet),  and  with  the  last  electron  occupying  ppa.  This 
leaves  the  ppa  bonding  state  half-occupied  and  able  to  ac¬ 
cept  one  more  electron.  Due  to  the  strong  bonding  nature  in 
the  CN  molecule,  this  ppa  state  has  a  low  energy.  When  the 
molecule  is  placed  in  GaAs,  the  ppa  state  lies  below  the 
VBM  and  always  becomes  fully  occupied.  At  the  same  time, 
the  insertion  of  CN  at  an  interstitial  site  also  creates  strain  in 
nearby  Ga-As  bonds,  leading  to  Ga-As  bond  extension  (or 
breaking).  The  broken  Ga-As  bonds  are  shown  in  dotted 
lines  in  Fig.  6.  CN;  is  stable  in  two  charge  states:  1+  and 
neutral.  We  initially  used  symmetric  configurations  in  our 
calculations  [Fig.  6(a)],  but  this  configuration  can  spontane¬ 
ously  relax  into  asymmetric  configurations  with  lower  ener¬ 
gies  [Fig.  6(b)  and  6(c)]  if  symmetry  breaking  is  allowed. 
The  energy  change  as  the  CN;  rotates  from  the  symmetric 
configuration  into  the  asymmetric  one  (for  1+  charge  state) 
calculated  by  using  the  NEB  method  is  shown  in  Fig.  7. 
Figure  7  shows  that  the  asymmetric  configuration  is  favored 


FIG.  6.  (Color  online)  Local 
atomic  geometry  of  interstitial  CN 
molecules  in  GaAs:  (a) 

CN^(sym),  (b)  CN^(asym),  and 
(c)  CN6(asym). 
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FIG.  7.  (Color  online)  Calculated  energy  of  CN-  +  as  the  CN 
molecule  rotates  from  a  symmetric  to  an  asymmetric  configuration. 
The  rotation  occurs  without  any  barrier  and  lowers  the  energy  by 
almost  0.6  eV. 

over  the  symmetric  one  by  —0.6  eV.  We  can  clearly  see  in 
the  asymmetric  configuration  [Fig.  6(c)]  that  only  one  of  the 
Ga-As  bonds  is  now  broken  (another  bond  is  extended  from 
a  typical  Ga-As  bond).  One  of  the  two  electrons  released 
from  the  broken  bond  goes  to  the  ppa  state  of  the  CN  mol¬ 
ecule,  whereas  another  one  is  removed,  resulting  in  the  1  + 
charge  state.  As  an  additional  electron  is  added  (increasing 
the  charge  state  of  the  defect  from  1+  to  neutral),  the  addi¬ 
tional  electron  goes  to  the  ppir*  level,  which  is  located  at 
approximately  0.5  eV  above  the  VBM  (see  Fig.  2).  This 
leads  to  the  +/0  transition  level  of  CN;  at  0.5  eV  above  the 
VBM.  The  occupation  of  the  ppir*  level  leads  to  the  increase 
in  the  CN  bond  to  1.25  A  and  is  reflected  in  an  enhanced 
binding  between  the  molecule  and  its  neighbors  [Fig.  6(c)] 
compared  to  the  1+  charge  state  [Fig.  6(b)}.  Because  the 
electron  density  of  the  ppir*  state  has  four  lobes  pointing 
outward,  which  is  similar  to  Fig.  3(b),  when  this  state  is 
occupied,  the  molecule  can  gain  energy  by  turning  the  lobes 
to  form  bonds  with  the  neighboring  atoms. 

Vibrational  frequencies  for  the  stretching  modes  of  CN;  in 
GaAs  and  GaP  are  shown  in  Tables  I  and  II.  The  frequency  is 
higher  for  the  CN;  in  the  1+  charge  state,  which  has  all  the 
bonding  states  occupied  and  the  antibonding  states  empty. 
The  frequency  of  1934  cm-1  (1897  cm-1  for  a  hypothetical 
symmetric  configuration)  is  comparable  to  that  of 
CN^(asym),  which  has  a  similar  electron  occupation  on  the 
molecule.  For  the  neutral  charge  state  CN;,  with  an  electron 
occupying  the  ppir*  state,  the  C-N  bond  length  increases  to 
1.25  A  and  the  vibrational  frequency  decreases  to 
1519  cm-1. 

C .  General  trends  of  the  bond  length  and  the  vibrational 
frequency 

Our  results  show  that  the  C-N  bond  distance  becomes 
larger  (and  the  vibrational  frequency  becomes  lower)  as 
more  electrons  are  added  into  the  ppir*  states.  Our  investi¬ 
gations  of  the  various  configurations  of  CN  molecules  in 
GaAs  and  GaP  allow  us  to  quantitatively  study  the  relation¬ 


FIG.  8.  (Color  online)  (a)  Calculated  vibrational  frequencies  and 
(b)  C-N  bond  distances  for  various  configurations  of  a  CN  molecule 
in  GaAs  and  GaP  as  a  function  of  charge  A q  inserted  into  the 
molecule,  which  are  measured  relative  to  a  neutral  free  CN  mol¬ 
ecule  with  nine  electrons.  The  solid  symbols  are  for  interstitial  CN 
and  the  open  symbols  are  for  substitutional  CN.  The  colors  code  the 
host  material  (red  for  GaAs  and  green  for  GaP),  whereas  the  sym¬ 
bols  are  used  to  distinguish  the  symmetric  (circle)  from  the  asym¬ 
metric  (triangle)  configurations.  The  CN  molecule  is  characterized 
as  triply,  doubly,  or  singly  bonded  for  Ag  =  -1,  -3,  and  -5,  respec¬ 
tively.  The  dashed  lines  are  linear  fits  to  all  of  the  data  points. 

ship  between  the  electron  occupation  (of  the  CN  molecule) 
and  the  bond  length  as  well  as  the  vibrational  frequency.  To 
do  this,  we  define  A q  as  the  additional  charge  being  added  to 
the  CN  molecule  relative  to  the  neutral  free  CN  molecule.  A 
free  CN  molecule  has  nine  electrons,  which  occupy  the  mo¬ 
lecular  states  up  to  half  of  the  ppa  state.  Since  this  state  is 
located  below  the  VBM  of  GaAs,  it  is  always  fully  occupied. 
This  means  that  at  least  one  electron  from  the  host  has  to  be 
added  to  the  molecule  (Ag  =  -1)  and  the  electron  occupation 
of  the  molecule  becomes  the  same  as  that  of  a  free  (triple- 
bonded)  CN“  ion.  For  the  CN  molecule  in  GaAs,  this  occu¬ 
pation  corresponds  to  the  lowest  charge  state,  i.e.,  CNA+  and 
CN-+.  As  the  Fermi  energy  is  raised,  electrons  are  inserted 
into  the  antibonding  states  {ppir*  states)  of  the  CN  molecule, 
thus  weakening  the  C-N  bond.  When  two  electrons  are  in¬ 
serted  into  the  antibonding  state  (i.e.,  CNAs),  making 
Ag  =  -3,  the  C-N  bond  is  reduced  to  a  double  bond  instead  of 
a  triple  bond.  When  four  electrons  are  inserted  into  the  anti¬ 
bonding  states  (i.e.,  CNA“),  making  Ag  =  -5,  the  C-N  bond  is 
reduced  to  a  single  bond.  Our  calculations  thus  show  that 
CNAs  can  have  its  bond  strength  varied  from  a  triple  bond 
down  to  a  single  bond  depending  on  the  electron  occupation. 
CN;,  on  the  other  hand,  can  only  exhibit  the  bond  strength  of 
a  triple  bond  (CN  - +)  or  a  triple  bond  with  one  electron  in  the 
antibonding  state  (CN^)  (which  can  be  considered  as  halfway 
between  a  double  bond  and  a  triple  bond). 

A  plot  of  the  C-N  bond  distance  as  a  function  of  A q  is 
shown  in  Fig.  8(b).  The  bond  distance,  especially  of  the  sym¬ 
metric  configurations,  clearly  show  a  linear  trend  with  A q. 
The  deviations  from  the  linear  trend  in  the  case  of  some  of 
the  asymmetric  configurations  can  be  attributed  to  the  inter¬ 
action  with  neighboring  host  atoms.  A  linear  fit  between  the 
bond  distance  and  A q  for  all  data  points  yields 
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dc_N  =  -0.0496A<?+ 1.143,  (3) 

where  dc_N  is  the  bond  distance  between  C  and  N  in  A.  The 
bond  length  increases  by  approximately  0.05  A  for  every 
electron  inserted. 

The  vibrational  frequency  decreases  as  the  bond  is  ex¬ 
tended.  The  relationship  between  the  frequency  and  A q  is 
shown  in  Fig.  8(a).  Again,  a  linear  trend  is  seen.  A  linear  fit 
yields 

(o  =  207 Ag  +  2057,  (4) 

where  a>  is  the  vibration  frequency  in  cm-1,  i.e.,  the  fre¬ 
quency  is  redshifted  by  approximately  200  cm-1  for  every 
electron  that  is  inserted.  Deviations  from  the  fit  (on  the  order 
of  100  cm-1)  occur  because  the  actual  frequency  also  de¬ 
pends  on  other  factors  such  as  the  detailed  local  geometry. 
Still,  the  linear  fit  in  Eq.  (4)  nicely  highlights  the  general 
trend,  even  if  it  does  not  provide  a  completely  accurate  pre¬ 
diction  of  the  frequency. 

D.  Discussion  and  comparison  with  experiment 

Ulrici  and  Clerjaud21  performed  a  polarized  infrared  ab¬ 
sorption  spectroscopy  on  C-doped  dilute-nitride  GaAs  and 
GaP.  They  observed  a  sharp  local  vibrational  mode  at 
2087.1  cm-1  (at  T=1  K)  in  GaP  and  a  similar  mode 
(2088.5  cm-1)  in  GaAs.  They  identified  it  as  a  CN  complex 
with  a  triple  bond  (between  C  and  N)  aligned  along  three 
equivalent  (100)  directions.  Their  C  and  N  assignments  were 
based  on  their  observation  of  additional  weak  peaks  at  203 1 
and  2048  cm-1,  which  exhibit  an  intensity  ratio  and  fre¬ 
quency  shift  that  match  the  natural  abundance  and  predicted 
frequency  shift  for  C13  and  N15  isotopes.  The  orientation  of 
the  molecule  was  confirmed  by  polarized  IR  absorption  un¬ 
der  different  types  of  uniaxial  stress. 

An  inspection  of  Tables  I  and  II  shows  that  the  measured 
frequency  of  2087  cm-1  is  most  consistent  with  our  calcu¬ 
lated  frequencies  for  triply  bonded  CN  molecules  with  A q 
=  -l.  A  more  detailed  identification  runs  into  some  compli¬ 
cations,  however.  There  are  four  configurations  with  triple 
bonds  covered  in  our  study:  CNA+(asym),  CNA+(sym), 

CN1+(asym),  and  CN1+(sym).  (1)  CNA+(asym)  provides  the 
best  matching  frequency  to  the  experiment  (lower  by 
36  cm-1,  which  is  well  within  the  computational  error  bar). 
However,  the  molecule  does  not  align  along  (100)  but,  in¬ 
stead,  prefers  to  align  in  a  direction  close  to  (110).  (2) 
CN^(sym)  has  the  correct  orientation;  however,  its  fre¬ 
quency  is  too  low  (by  309  cm-1),  and,  even  more  impor¬ 
tantly,  it  is  not  energetically  stable  and  spontaneously  rotates 
to  CNA+(asym),  lowering  the  energy  by  0.27  eV.  Therefore, 
we  can  rule  out  CNA+(sym).  (3)  CN1+(asym)  is  a  possible 
candidate.  Its  calculated  vibrational  frequency  is  not  too  far 
from  the  observed  value  (lower  by  154  cm-1).  However,  its 
orientation  also  deviated  from  (100)  by  16°.  In  addition,  the 
interstitial  configurations  are  energetically  unfavorable  com¬ 
pared  to  those  of  CNA+(asym).  (4)  CN1+(sym)  is  an  unlikely 
candidate  because  it  is  not  energetically  stable.  It  can  spon¬ 
taneously  rotate  to  the  asymmetric  configuration  and  gain  0.6 
eV.  In  addition,  its  vibrational  frequency  is  even  lower  than 


that  of  the  asymmetric  configuration,  widening  the  disagree¬ 
ment  with  experiment  to  191  cm-1. 

Our  first-principles  calculations  of  vibrational  frequencies 
for  CN  molecules  in  various  configurations  in  GaAs  and  GaP 
therefore  generally  support  Ulrici  and  Clerjaud’ s21  assign¬ 
ment  of  the  2087  cm-1  mode  to  a  triply  bonded  CN  mol¬ 
ecule.  However,  our  calculations  also  show  that  the  triply 
bonded  CN  molecule  (in  either  the  substitutional  or  intersti¬ 
tial  configurations)  does  not  symmetrically  align  in  a  tetra¬ 
hedral  site  surrounded  by  four  Ga  atoms.  As  a  result,  neither 
CNAs  nor  CN;  have  the  C-N  axis  aligned  in  the  (100)  direc¬ 
tion,  as  proposed  by  Ulrici  and  Clerjaud.21  However,  given  a 
rather  low  rotation  barrier  of  only  0.27  eV  (Fig.  4),  it  is 
possible  that  CNA+(asym)  is  constantly  rotating,  resulting  in 
an  average  orientation  in  the  (100)  direction  observed  in  the 
experiment. 

IV.  CONCLUSIONS 

We  have  presented  first-principles  results  for  CN  mol¬ 
ecules  in  GaAs  and  GaP.  The  study  covered  two  possible 
lattice  locations:  (1)  CN  molecules  substituting  for  anions 
(CNAs(orP))  and  (2)  interstitial  CN  molecules  at  the  Td  site 
surrounded  by  four  Ga  atoms  (CN;).  All  possible  charge 
states  and  various  orientations  were  considered.  The  calcu¬ 
lated  formation  energies  show  that  the  molecule  favors  sub¬ 
stituting  for  the  anion  site  over  the  interstitial  configuration 
for  all  equilibrium  growth  conditions  with  a  margin  of  at 
least  1  eV.  The  calculations  predict  the  CN  molecule  to  pro¬ 
duce  a  level  with  a  strong  ppir*  molecular  orbital  character¬ 
istic  at  approximately  0.5-0. 8  eV  above  the  VBM.  In  p-type 
conditions,  where  the  Fermi  level  is  located  below  this  ppir* 
level,  the  molecules  are  triply  bonded  and  form  donor  defect 
centers  (double  donor  CNA+  and  single  donor  CN1+).  The 
calculations  with  full  relaxation  show  that  triply  bonded 
CNA+  and  CN1+  do  not  symmetrically  orient  in  the  Ga  tetra¬ 
hedron  but,  instead,  are  tilted  in  order  to  gain  better  interac¬ 
tions  with  neighboring  Ga  atoms.  The  calculated  vibration 
frequencies  of  the  triply  boned  CNA*  and  CN1+  are  in  rea¬ 
sonable  agreement  with  the  measurement  by  Ulrici  and  Cler¬ 
jaud  supporting  their  identification  that  the  bonding  is  triple¬ 
bond  type.  Although  neither  CNAs  nor  CN;  have  the  C-N 
bond  oriented  along  the  (100)  directions,  as  proposed  by 
Ulrici  and  Clerjaud,21  the  low  rotation  barrier  of  CNA+  makes 
it  possible  that  the  molecule  might  be  constantly  rotating, 
leading  to  an  average  orientation  in  the  (100)  direction.  At 
higher  Fermi  levels,  the  CN1+  can  accept  an  electron  and 
becomes  CN 9,  whereas  CNA+  can  accept  one,  two,  three,  or 
four  electrons  (depending  on  the  Fermi  level)  and  becomes 
CNA+,  CNAs  (metastable),  CNA“,  or  CNA“,  respectively.  We 
found  that  the  C-N  bond  length  and  vibrational  frequency 
linearly  change  (to  a  good  approximation)  with  the  number 
of  electrons  added  into  the  antibonding  states  of  the  mol¬ 
ecule. 
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Current  controversy  on  the  binding  sites  of  H+  in  ZnO  can  be  explained  by  first-principles  calculations. 

Previous  infrared  measurements  from  different  groups  indicate  different  H  sites,  either  at  the  bond  center  (BC) 
site  with  a  stretch  frequency  to=3611  cm-1  or  antibonding  (AB)  site  with  to=3326  cm-1.  This  was  puzzling 
because  the  BC  site  has  lower  energy  by  0.2  eV.  Here,  we  show  that  calcium,  isovalent  to  Zn  and  found  only 
in  samples  with  the  3326  cm-1  mode,  binds  H  at  the  AB  site.  Large  spatial  undulation  of  charge  explains  the 
unexpected  large  binding  between  isovalent  Ca  and  charged  H+  of  0.7  eV. 

DOI:  10. 1103/PhysRevB.78. 113203  PACS  number(s):  61.72.Bb,  63.20.Pw 


Hydrogen  in  ZnO  has  attracted  a  lot  of  attention  in  the  last 
decade.  In  2000,  Van  de  Walle  proposed  that  H  prefers  to 
form  a  strong  O-H  bond  with  lattice  O  and  acts  exclusively 
as  a  donor,  H+.1  This  behavior  is  different  from  H  in  almost 
all  other  semiconductors  in  which  H  is  amphoteric,  counter¬ 
acting  the  dominant  dopant.2  ZnO  is  a  potential  superior 
electronic  material  for  next-generation  blue  optoelectronic 
applications.  However,  due  to  the  lacking  of  high  quality 
/?-type  samples,  fabrication  of  ZnO-based  optoelectronic  de¬ 
vices  is  still  not  possible.  The  prediction  that  H  is  an  exclu¬ 
sive  donor  stirred  up  considerable  research  activities  because 
it  may  explain  the  difficulty  in  p  doping.3-11 

Binding  site  is  one  of  the  most  fundamental  properties  of 
any  impurity.  Although  it  is  clear  that  H+  prefers  tc>o  bind 
strongly  with  O  at  a  bond  length  of  approximately  1  A  and 
with  a  vibration  frequency  around  3500  cm-1,  the  details  is 
completely  controversial:  McCluskey  et  al.6  observed  an  IR 
peak  at  3326  cm-1  in  samples  provided  by  Cermet,  while 
Lavrov  et  al.12  observed  a  different  peak  at  3611  cm-1  in 
samples  by  Eagle-Picher.  Polarized  IR  showed  that  the 
3326  cm-1  peak  is  associated  with  an  oscillator  orientated  at 
112°  angle  with  respect  to  the  c  axis,  whereas  those  that  give 
the  3611  cm-1  peak  are  parallel  to  the  axis.  Later  study  by 
Shi  et  al.1  on  both  samples  showed  that  the  Cermet  samples 
have  a  strong  3326  cm-1  peak  and  a  weak  3611  cm-1  peak 
but  the  Eagle-Picher  samples  do  exactly  the  opposite.  This 
indicates  that  there  exist  two  forms  of  O-H  in  these  two  types 
of  samples.  Because  the  splitting  3611-3326=285  cm-1  is 
comparable  to  the  calculated  splitting  between  the  BC  and 
AB  site  adjacent  to  O  (AB0),  the  3611  cm-1  peak  has  been 
assigned  to  BCy  and  the  3326  cm-1  peak  to  AB0i,13  where  II 
and  J_  denote  the  orientation  of  the  O-H  bond  with  respect  to 
the  c  axis  [see  Fig.  1(a)]. 

However,  independent  calculations  by  different 
groups12,13  showed  consistently  that  the  BC||  site  has  lower 
formation  energy  than  any  other  site  by  0.2  eV,  which,  in 
comparison  with  the  experimental  temperature  of  4  K,  is 
quite  large.  More  importantly,  it  is  difficult  to  understand 
why  the  AB0i  configuration  can  dominate  in  some  samples 


but  not  in  others.  An  alternative  explanation  is  that  one  (or 
both)  of  the  signals  may  come  from  an  H  complex  rather 
than  just  an  isolated  H+.  Recently,  McCluskey  and  Jokela14 
studied  the  trace-amount  elements  in  both  samples  by  using 
delayed  gamma  neutron  activation  analysis  (DGNAA)  and 
secondary-ion  mass  spectrometry  (SIMS).  They  found  that 
the  Cermet  samples,  grown  by  a  pressurized  melt-growth 


(a) 


FIG.  1.  (Color  online)  [(a)  and  (b)]  Schematics  in  the  (11-20) 
plane  of  the  four  H+  sites  in  Ca-free  ZnO  and  near  a  Ca  impurity. 
Calculated  atomic  relaxations  of  Ca-AB0i  (c)  in  the  absence  and 
(d)  in  the  presence  of  H+.  Oxygen  atoms  adjacent  to  the  Ca  are 
marked  by  darker  blue.  Dashed  circles  and  dashed  lines  indicate  the 
original  bulk  positions  and  bonds.  Distances  are  given  either  in  unit 
of  angstrom  or  in  percentage  with  respect  to  that  of  bulk  ZnO. 
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TABLE  I.  Calculated  formation  energy,  O-H  bond  length,  and  stretch  frequencies  for  H+  in  ZnO  without 
(upper  panel)  and  with  (lower  panel)  Ca.  The  energy  zero  is  that  of  Ca-free  BCy  site.  Highlighted  rows  are  the 
ground- state  configurations. 


Site 

A£(eV) 

^O-H 

A 

Vibration  frequency, 

w(cm  ') 

Present 

Ref.  13 

Present 

Ref.  13 

Present 

Ref.  13 

Measured 

BCy 

OlOO 

OlOO 

0l985 

Ol906 

3421 

3475 

3611 

BC± 

0.14 

0.15 

0.985 

0.982 

3505 

3519 

ABqu 

0.19 

0.17 

1.001 

1.003 

3097 

3116 

AB0i 

0.15 

0.14 

1.004 

1.003 

3109 

3154 

Ca-BCy 

0.54 

0.967 

3934 

Ca-BC± 

0.38 

0.979 

3484 

Ca-AB0|| 

-0.29 

0.995 

3236 

Ca-AB0± 

-0l49 

0l993 

3207 

3326 

Reference  12. 
bReference  6. 


process,15  contain  isovalent  Ca  impurity,  whereas  in  the 
Eagle-Picher  samples,  which  were  grown  by  chemical  vapor 
transport,16  the  amount  of  Ca  impurity  is  below  the  detection 
threshold. 

In  this  Brief  Report,  we  show  that  Ca  may  play  a  pivotal 
role  in  switching  the  H  from  the  BCy  site  to  the  AB0j_  site. 
Our  calculation  shows  that,  with  Ca,  the  AB0i  site  is  0.5  eV 
more  stable  than  the  ground-state  BCy  site  in  pure,  Ca-free 
ZnO.  Our  results  are  consistent  with  the  experimental  finding 
that,  for  Cermet  samples  with  Ca,  the  3326  cm-1  line  is 
observed  with  the  transition  moment  aligned  at  112°  from 
the  c  axis,  whereas  for  Eagle-Picher  samples  without  Ca,  the 
3611  cm-1  line  with  a  transition  moment  aligned  along  the  c 
axis  is  observed.  The  fact  that  an  isovalent  impurity  can  bind 
H+  so  strongly  by  0.5  +  0.2  =  0.7  eV  is  unexpected  and  hence 
surprising.  Our  analysis  shows  that  this  could  be  a  general 
property  of  ionic  semiconductors  and  solids  for  which  the 
electrostatic  monopole  (such  as  H+)  and  dipole  (due  to  spa¬ 
tial  undulated  local  field)  interactions  can  be  comparable  to 
the  monopole-monopole  interactions  (such  as  within  a 
donor- acceptor  pair)  in  more  covalent  materials. 

We  used  the  density-functional  theory17  within  the  local- 
density  approximation  (LDA)  and  the  projector  augmented 
wave  method18,19  as  implemented  in  the  VASP  code.20,21  Zinc 
3d  states  are  treated  as  valence  states.  The  cut-off  energy  for 
the  plane-wave  expansion  is  400  eV.  We  used  a  supercell 
approach  (96-atom  wurtzite  cell)  with  the  Monkhorst-Pack 
&-point  mesh  for  the  Brillouin-zone  integration  (2X2X2). 
For  charged  defect,  we  used  the  jellium  background  approxi¬ 
mation.  All  the  atoms  are  allowed  to  relax  to  their  equilib¬ 
rium  positions.  To  calculate  the  O-H  frequencies,  we  em¬ 
ployed  the  approach  in  Ref.  13  to  include  anharmonic 
effects.  For  free  H20  molecule,  this  yields  a  symmetric 
stretch  frequency  of  3540  cm-1,  which  is  lower  than 
experiment22  by  coer=117  cm-1  (an  error  of  —3%).  To  cor¬ 
rect  such  a  systematic  error,  ojer  has  been  added  to  all  the 
calculated  results  before  comparing  with  experiments. 

Table  I  shows  the  calculated  formation  energies  A E  [rela¬ 
tive  to  the  ground-state  (BCy)  energy  in  Ca-free  ZnO],  O-H 
bond  lengths  d0_ H,  and  stretch  frequencies  co  for  H+  in  ZnO 


without  Ca  (upper  panel)  and  with  Ca  (lower  panel).  To  be 
consistent,  the  upper  panel  are  the  results  of  current  96-atom 
cell  calculations,  which  are  qualitatively  the  same  as  previ¬ 
ous  calculations,  e.g.,  d0_ H  and  co  agree  with  our  previous 
results  to  within  1%  and  3%,  respectively.13  With  Ca,  we 
denote  the  complexes  as  Ca-BC||,  Ca-BC±,  Ca-AB0||,  and 
Ca-ABoi,  respectively,  with  their  atomic  structures  sche¬ 
matically  shown  in  Fig.  1(b).  The  formation  energies  of  the 
Ca-containing  complexes  are  defined  as 

=  £tot(Ca-H+)  +  £tot(bulk)  -  EtjCaZn)  -  £tot(BC||), 

(1) 

where  Etot(a)  is  the  total  energy  of  a  supercell  containing  the 
complex  (or  impurity)  a.  Physically,  A E  is  the  energy  gain  or 
loss  by  bringing  a  H+  from  the  lo west-energy  BCy  site  in 
Ca-free  ZnO  to  the  respective  sites  near  the  Ca  in  Fig.  1(b). 
A  negative  A E  implies  that  the  H+  is  bounded,  whereas  a 
positive  A E  implies  that  the  H+  is  unbounded,  to  the  Ca. 
Even  when  A E  is  positive,  the  site  is  still  metastable  due  to 
the  diffusion  barriers  of  H+  in  ZnO. 

As  Ca  is  isovalent  to  Zn,  on  a  first  look,  it  is  puzzling  why 
a  charged  H+  would  like  to  bind  to  a  neutral  Ca  with  a  rather 
sizable  energy — a  behavior  usually  found  either  in  oppo¬ 
sitely  charged  donor-acceptor  pairs  or  in  charge  neutral  but 
large-strain  compensated  defect  pairs.  Our  calculation  shows 
that  both  substitutional  Ca  and  interstitial  H+  cause  outward 
relaxations  to  the  surrounding  Zn  and  O  atoms.  Therefore, 
neither  of  the  above  can  be  the  cause  for  the  unusually  large 
binding  energy.  It  is  even  more  puzzling  why  the  site  pref¬ 
erence  is  reversed  from  the  BC  site  to  the  AB0  sites. 

To  understand  the  puzzles,  we  note  first  that  H+  in  ZnO 
binds  to  the  oxygen  atoms  rather  than  any  of  the  cation  at¬ 
oms.  Second,  ZnO  is  a  highly  ionic  semiconductor  with  an 
ionicity  of  0.616  in  Phillips’  scale.23  In  other  words,  although 
ideal  ZnO  is  a  homogeneous  bulk  material,  there  is  a  large 
spatial  undulation  of  the  electric  field  at  the  atomic  level. 
Therefore,  the  defect  and  impurity  physics  that  is  already 
well  established  for  covalent  semiconductors  such  as  Si  or 
GaAs  may  not  necessarily  apply  here.  Rather,  it  is  how  much 
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FIG.  2.  (Color  online)  Radial  distributions  of  neighboring  atoms 
surrounding  an  H+  at  different  sites:  (a)  BCy,  (b)  Ca-BCy,  (c)  AB0^, 
and  (d)  Ca-AB0i.  Solid  lines  highlight  the  local  bonds:  O-H-Zn  (or 
O-H-Ca)  for  the  BC  sites  and  H-O-Zn  (or  H-O-Ca)  for  the  AB0 
sites.  Dashed  circles  indicate  the  radial  distances  from  the  H+  cen¬ 
ter.  For  simplicity,  atoms  in  the  shell  are  evenly  distributed.  In  (b) 
and  (d),  shells  with  enhanced  charge  with  respect  to  (a)  and  (c)  are 
highlighted  by  halo:  blue  for  enhanced  negative  and  red  for  en¬ 
hanced  positive. 

the  impurity  affects  the  undulation  of  the  field  that  holds  the 
key  to  defect  and  impurity  physics  in  ionic  semiconductors. 
The  electronegativity  of  Ca,  Zn,  and  O,  on  the  Pauling  scale, 
is  *Ca=1.00,  /YZn=  1  -65,  and  ^°=3.44.  Therefore,  Ca  has  a 
stronger  tendency  to  donate  its  electrons  to  neighboring  O 
atoms  than  Zn  does.  This  stronger  charge  transfer  makes  the 
O  atoms  adjacent  to  the  Ca  atom  stronger  negative  centers 
relative  to  other  O  atoms.  As  such,  qualitatively  speaking, 
the  H+  should  bind  to  these  O  atoms  more  strongly.  Also, 
because  Ca  is  a  much  stronger  positive  charge  center  than 
Zn,  a  H+  would  avoid  the  BC  sites  adjacent  to  Ca,  making 
the  adjacent  AB0  sites  energetically  more  favorable. 

To  further  the  discussion,  next  we  apply  a  Bader 
analysis.24,25  Our  calculations  show  that  in  bulk  ZnO,  each 
Zn  donates  1.13  electrons  to  O.  In  other  words,  Zn  can  be 
considered  as  a  positive  center  with  q=  +  I.13e  whereas  O 
can  be  considered  as  a  negative  center  with  q=- l.\3e.  This 
explains  why  an  H+  prefers  to  stay  close  to  the  O,  instead  of 
Zn.  Bader  analysis  of  CaZn  further  shows  that  Ca  donates 
lA2e,  which  is  0.29e  more  than  what  Zn  does.  The  extra 
donation  is  almost  equally  distributed  among  the  four  O 
nearest  neighbors:  0.06^  to  the  O  in  the  [0001]  direction,  and 
0.01  e  to  each  of  the  remaining  three  O. 

To  understand  the  switching  in  site  preference  from  BC  in 
bulk  ZnO  to  AB0  near  CaZn,  we  need  to  examine  next- 
nearest  neighbors  (~2  A)  where  the  difference  between  the 
BC  and  AB0  sites  can  be  readily  seen.  Figure  2  shows,  in  a 
schematic  shell  model,  the  radial  distributions  of  neighboring 
atoms  of  H.  It  so  happens  that  before  introducing  the  Ca,  the 
four  sites  in  the  upper  panel  of  Table  I,  BCy,  BC±,  AB0||,  and 
AB0i,  have  reasonably  close  formation  energies  to  within 


0.2  eV.  This  enables  us  to  interpret  the  changes  upon  Ca 
substitution  primarily  in  terms  of  a  redistribution  of  the 
charge.  Because  the  BC  sites  are  closer  to  the  Ca  [e.g.,  Fig. 
2(b)],  they  become  energetically  unfavorable.  In  contrast,  the 
AB0  sites  are  much  further  away  from  the  Ca  [see,  Fig. 
2(d)].  As  a  result,  the  AB0  sites  benefit  more  from  the  charge 
redistribution,  forming  a  stronger  0-H+  bond  yet  staying 
away  from  the  positive  Ca  so  as  not  to  be  strongly  affected 
by  it.  Moreover,  because  H  in  the  BC  sites  cuts  the  Ca-0 
bonds,  it  effectively  blocks  charge  transfer  from  Ca  to  O.  As 
a  matter  of  fact,  Bader  charge  on  the  first  shell  O  atom  in 
Fig.  2(b)  is  the  same  as  the  Ca-free  case.  Hence,  not  only 
does  the  Ca-BC  sites  become  less  favorable  than  the  Ca-AB0 
sites,  the  formers  are  positive  in  energy  with  respect  to  the 
BCy  site  away  from  the  Ca.  In  other  words,  H+  will  not  bind 
to  the  Ca  at  the  Ca-BC  sites. 

Next,  we  estimate  the  equilibrium  distribution  of  H+ 
among  the  various  lattice  sites  a.  By  constructing  the  parti¬ 
tion  functions,  we  can  derive  the  equations  relating  the  con¬ 
centration  of  H+(a),  [H+(a)],  to  the  overall  chemical  poten¬ 
tial  of  H+,  /xH+.26  Because  we  consider  only  H+,  Fermi-level 
dependence  of  the  various  sites  a  is  the  same.  Vibrational 
energy  contributions  for  these  sites  are  similar,  as  well.  Their 
contributions  can  thus  be  expressed  as  E0  and 

fa.  -  fa.  —  fa  —  Mia) +  kT  ln{  j(  } , 

(2) 

where  A E(a)  is  the  formation  energy  of  H+(cr)  in 
Table  I,  n(a)  is  the  density  of  available  lattice  sites  for  a 
specific  H+(a),  and  j  is  the  number  of  symmetry  degeneracy. 
For  example,  for  H+  at  BCy  site,  we  have 

/zH+=0.0  eV+krin{[H+(BCy)]/[0]-[H+(BCy)]},  and  for 
H+  at  Ca-ABoi  site,  we  have  /IH+=-0.49  eV 

+kT  ln([H+(Ca-AB0i)]/3{[Ca]-[H+(Ca-AB0x)]}).  To 
solve  the  set  of  equations,  we  apply  the  equilibrium  con¬ 
straint  that  /xH+  is  the  same  for  all  the  H+(a).  We  also  require 
that  the  total  concentration  of  H+  satisfies  [H+]  =  2[H+(a)], 
whose  value  is  an  input  from  experimental  conditions. 

For  Cermet  samples,  the  concentration  of  Ca  is  measured 
to  be  [Ca]  =  4  X  1016  cm-3  (Ref.  14)  and  the  total  concentra¬ 
tion  of  H+  is  expected  to  be  at  most  in  the  same  order  of 
magnitude  (so  we  use  [H+]  =  4  X  1016  for  the  illustration  pur¬ 
pose).  Figure  3(a)  shows  the  normalized  equilibrium  concen¬ 
trations  as  a  function  of  temperature.  Concentrations  with 
less  than  5%  of  total  [H+]  are  ignored  in  the  figure.  We  see 
that  the  majority  (over  75%)  of  H+  is  on  the  Ca-AB0i  sites 
for  temperatures  below  400  K  (or  127  °C).  It  has  been  pre¬ 
dicted  that  H+  diffuses  readily  at  room  temperature  due  to  its 
low  migration  barrier  of  less  than  0.5  eV  (Refs.  3  and  27)  to 
equilibrate  the  H+  among  the  various  lattice  sites.  For  Eagle- 
Picher  samples,  the  concentration  of  Ca  is  outside  experi¬ 
mental  detection  limit  and  is  hence  much  less  than  that  of 
H+.  Figure  3(b)  shows  the  results  for  such  a  case.  Here  as 
expected,  the  majority  (over  75%)  of  H+  is  on  the  BCy  sites 
for  temperatures  less  than  600  K. 

The  calculated  vibrational  signatures  of  the  H+  sites  are 
also  in  agreement  with  experiments.  For  example,  the  stretch 
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FIG.  3.  (Color  online)  Normalized  equilibrium  concentrations  of 
[H+(<T)]  at  the  various  sites  a  as  a  function  of  temperature,  (a)  is  the 
case  when  [Ca]  =  4  X  1016  cm-3  and  (b)  is  the  case  when  there  is  no 
Ca.  A  total  [H+]=4X  1016  cm-3  was  assumed  in  the  calculations. 

frequencies  of  3421  cm-1  for  BCy  and  3207  cm-1  for 
Ca-AB0i  agree,  within  computational  uncertainty,  with  the 
experimental  lines  at  3611  (Eagle-Picher)  and  3326  cm-1 
(Cermet).  Note  that  although  the  absolute  values  here  differ 
from  experiments  by  about  100-200  cm-1,  the  difference 


between  the  two  experimental  modes  of  3611-3326 
=  285  cm-1  is  in  a  good  agreement  with  that  between  the  BC 
and  AB0  configurations  of  3421-3207  =  214  cm-1.  The  BQ 
configuration  has  the  oscillator  oriented  along  the  c  axis, 
which  is  consistent  with  the  experimental  observation  for  the 
Eagle-Picher  samples.  The  Ca-AB0i  configuration  has  the 
oscillator  orientated  at  107°  with  respect  to  the  c  axis,  which 
is  consistent  with  the  polarized  IR  measurement  for  the  Cer¬ 
met  samples  with  transition  moment  aligned  at  112°. 

In  conclusion,  the  large  atomic-level  spatial  undulation  of 
electric  field  makes  the  defect  and  impurity  physics  in  ionic 
materials,  in  particular  in  ZnO,  qualitatively  different  from 
that  of  covalent  semiconductors.  Larger  than  expected  bind¬ 
ing  energy  of  0.7  eV  can  result  between  charge  neutral,  is- 
o valent  Ca,  and  H+  impurities.  This  tilts  the  energy  balance 
between  the  BC  and  AB0  sites  to  favor  AB0i  in  the  pres¬ 
ence  of  Ca.  The  calculated  vibrational  signatures  of  H+  are 
also  consistent  with  IR  measurements.  This  may  explain  the 
current  puzzle  about  the  H+  sites  in  ZnO. 
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First-principles  calculations  show  unexpectedly  large  binding  (0.7  eV)  between  positive  H+  at 
oxygen  antibonding  site  (A^q)  and  Ca,  which  replaces  iso  valent  Zn  in  ZnO.  It  is  explained  by  a 
large  spatial  undulation  of  the  host  charge  and  the  change  in  the  charge  due  to  Ca  replacement  of  the 
Zn.  Our  findings  are  consistent  with  a  number  of  infrared  measurements.  ©  2009  American  Vacuum 
Society.  [DOI:  10.1116/1.3119679] 


I.  INTRODUCTION 

ZnO  is  a  potential  electronic  material  for  next-generation 
blue  optoelectronic  applications.  However,  due  to  the  lack  of 
good  quality  p-type  materials,  fabrication  of  ZnO-based  op¬ 
toelectronic  devices  is  still  not  possible.  Much  effort  has 
been  made  both  experimentally  and  theoretically  trying  to 
explain  the  p-type  doping  difficulty.  Hydrogen  may  be  a 
killer  as  first  proposed  by  Van  de  Walle1  in  2000.  Hydrogen 
is  a  very  common  impurity  in  semiconductors  because  it  is 
abundant  in  air  and  is  widely  used  in  various  growth  tech¬ 
niques.  Hydrogen  is  also  amphoteric:  in  p-type  semiconduc¬ 
tors,  it  acts  as  a  donor,  but  in  ft-type  semiconductors,  it  acts 
as  an  acceptor.  In  other  words,  it  is  always  opposite  to  the 
prevailing  conductivity  of  the  host.  However,  hydrogen  is 
very  special  in  ZnO.  It  prefers  to  form  a  strong  O-H  bond 
with  lattice  O  and  acts  exclusively  as  a  donor,  H+.1,2  There¬ 
fore,  H  in  ZnO  has  attracted  a  lot  of  attention.3-11 

Despite  the  initial  suggestion,  the  most  stable  configura¬ 
tion  of  hydrogen  in  ZnO  is  still  controversial.  Independent 
predictions  by  different  groups12,13  showed  consistently  that 
bond  center  BC\\  site  [see  Fig.  1(a)]  is  the  most  stable  con¬ 
figuration  with  formation  energy  less  than  any  other  sites  by 
at  least  0.2  eV.  Here,  we  use  II  and  _L  to  denote  the  orienta¬ 
tion  of  the  O-H  bond  with  respect  to  the  c-axis  of  ZnO. 
However,  experimental  observations  deviate  from  each 
other:  McCluskey  et  al.6  observed  an  IR  peak  at  3326  cm-1 
in  samples  provided  by  Cermet,  while  Lavrov  et  al.12  ob¬ 
served  a  different  peak  at  3611  cm-1  in  samples  by  EaglePi- 
cher.  Polarized  IR  showed  that  the  3326  cm-1  peak  is  asso¬ 
ciated  with  an  oscillator  oriented  at  112°  angle  with  respect 
to  the  c-axis,  whereas  the  3611  cm-1  peak  is  parallel  to  the 
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c-axis.  Later  study  by  Shi  et  al.1  on  both  samples  showed 
that  the  Cermet  samples  have  a  strong  3326  cm-1  peak  and  a 
weak  3611  cm-1  peak,  but  the  EaglePicher  samples  have 
exactly  the  opposite.  This  indicates  that  there  exist  two  forms 
of  O-H  in  these  samples.  Furthermore,  the  splitting 
3611-3326=285  cm-1  is  comparable  to  the  calculated  split¬ 
ting  between  the  BC  site  and  the  AB  site  adjacent  to  oxygen 
(A50),  the  3611  cm-1  peak  has  been  assigned  to  BC\\ , 
whereas  the  3326  cm-1  peak  to  AB01.13  The  dominance  of 
the  ABQ1  site  in  the  Cermet  samples  contradicts  with  theory 
as  a  0.2  eV  difference  is  considerably  larger  than  the  thermal 
energy  at  the  experimental  temperature  of  4  K. 

Recently,  McCluskey  and  Jokela,14  by  using  delayed 
gamma  neutron  activation  analysis  and  secondary  ion  mass 
spectrometry,  found  that  the  Cermet  samples  contain  consid¬ 
erable  amount  of  Ca,  whereas  the  EaglePicher  samples  do 
not.  This  raises  an  important  question  as  to  whether  the  un¬ 
intentional  impurity,  although  iso  valent  to  zinc,  can  cause 
size  switch  in  some  unexpected  way. 

In  this  study,  we  show  that  Ca  has  strong  Coulomb  bind¬ 
ing  with  H+  in  the  AB0 ±  site  despite  the  fact  that  Ca  is 
iso  valent  to  Zn.  Our  finding  is  consistent  with  hydrogen  in 
the  Cermet  samples,  i.e.,  the  3326  cm-1  IR  peak  dominates 
with  the  oscillator  aligned  112°  from  the  c-axis.  A  new  bind¬ 
ing  mechanism  based  on  charge  transfer  is  proposed  to  ex¬ 
plain  the  unexpected  strong  binding.15 

II.  METHOD 

Our  calculations  use  first-principles  density  functional 
theory16  within  the  local  density  approximation  and  the  pro- 
jector  augmented  wave  method,  ’  as  implemented  in  the 
VASP  code.19,20  Zn  3d  states  are  treated  as  the  valence  states. 
The  cutoff  energy  for  the  plane  wave  expansion  is  400  eV. 


1601  J .  Vac.  Sci.  Technol.  B  27(3),  May/J  un  2009  1071-1023/2009/27(3)/1601/3/$25.00  ©2009  American  Vacuum  Society  1601 
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Fig.  1.  (Color  online)  (a)  and  (b)  are  schematics  of  the  four  H+  sites  in  a 
(11-20)  plane  in  the  Ca-free  ZnO  and  near  a  Ca  impurity,  (c)  and  (d)  are  the 
calculated  atomic  relaxations  for  Ca  and  Ca-H+  at  the  Ca -AB0±  site,  re¬ 
spectively.  Oxygen  atoms  adjacent  to  the  Ca  are  marked  by  darker  blue. 


The  calculated  equilibrium  lattice  constant  of  ZnO  is  only 
1 .4%  shorter  than  experiment,  and  the  calculated  heat  of  for¬ 
mation  of  3.51  eV  agrees  with  the  experimental  value  of  3.60 
eV.  We  use  48-  and  96-atom  wurtzite  supercells  to  study  the 
Ca-H  pair  with  the  Monkhorst-Pack  &-point  mesh  for  the 
Brillouin  zone  integration  (3X3X2  for  48-atom  cell  and  2 
X2X2  for  96-atom  cell).  Unless  noted  otherwise,  the  nu¬ 
merical  values  presented  are  for  the  96-atom  cell  calcula¬ 
tions.  For  charged  defect,  we  use  the  jellium  background 
approximation.  All  the  atoms  are  relaxed  to  their  equilibrium 


positions.  To  calculate  the  O-H  frequencies,  we  employ  the 
approach  in  Ref.  13  to  include  anharmonic  effects.  For  free 
H20  molecule,  this  yields  a  symmetric  stretch  frequency  of 
3540  cm-1,  which  is  lower  than  experiment21  by  o)ER 
=  117  cm-1  (an  error  of  —3%).  To  correct  for  such  a  system¬ 
atic  error,  we  add  &>ER  to  all  the  calculated  frequencies  before 
comparing  them  with  experiments. 


III.  RESULTS  AND  DISCUSSION 

Figure  1  shows  the  H-related  atomic  configurations  with 
and  without  the  Ca.  In  agreement  with  previous  calculations, 
we  also  find  that  the  BC\\  site  is  the  most  stable  configuration 
in  Ca-free  ZnO.  Hence,  we  take  the  BC\\  site  as  the  energy 
reference  to  define  our  formation  energies, 

=  £tot(Ca-H+)  +  Etot(bulk)  -  £tot(CaZn)  -  Etot(BQ , 

(1) 

where  Et0 t(a)  is  the  total  energy  of  a  supercell  containing  the 
impurity  a.  In  a  way,  A E  is  the  energy  required  to  transfer  a 
BC\\  site  H+  from  Ca-free  ZnO  to  a  Ca-related  site  in  Ca- 
containing  ZnO.  A  negative  A E  implies  that  H+  prefers  to  go 
to  Ca,  whereas  a  positive  A E  implies  that  H+  does  not  prefer 
to  go,  or  in  other  words,  the  Ca-H  pair  is  only  metastable. 

Table  I  shows  the  calculated  formation  energy  and  the 
O-H  distance  for  Ca-free  and  Ca-containing  ZnO.  From  the 
table,  we  can  see  that  (1)  the  AB01  site  next  to  a  Ca  is  the 
most  stable  configuration.  This  is  consistent  with  experi¬ 
ments  on  the  Cermet  samples.  (2)  Surprisingly,  at  the 
Csl-ABq1  site,  the  Ca  atom  and  H+  ion  have  a  strong  binding 
(0.5 +  0.2= 0.7  eV,  which  is  almost  comparable  to  the  oppo¬ 
sitely  charged  donor  and  acceptor  pairs).  Our  calculation 
shows  that  both  substitutional  Ca  and  interstitial  H+  cause 
outward  relaxations  to  the  surrounding  Zn  and  O  atoms. 
Thus,  the  strong  binding  here  cannot  be  strain  related.  If  it  is 
neither  Coulomb  (for  opposite  charge)  nor  lattice  strain,  we 
need  a  new  mechanism  to  explain  our  results. 

First  of  all,  we  note  that,  based  on  the  Pauling  scale,  the 
electronegativities  of  Ca,  Zn,  and  O  are  ^Ca=1.00, 

=  1.65,  and  ^°  =  3.44,  respectively,  and  ZnO  is  an  ionic  semi- 


Table  I.  Calculated  hydrogen  formation  energy  (relative  to  H  on  the  BCl{  site  in  the  Ca-free  ZnO)  and  O-H 
bond  lengths.  Highlighted  rows  are  for  the  respective  ground-state  configurations  without  and  with  Ca. 


A  E 
(eV) 

^O-H 

(A) 

Site 

96-atom  supercell 

48-atom  supercell 

96-atom  supercell 

48-atom  supercell 

BC\\ 

0  00 

0  00 

0  985 

0  986 

BCl 

0.14 

0.16 

0.985 

0.988 

A£o|| 

0.19 

0.18 

1.001 

1.010 

ABOL 

0.15 

0.17 

1.004 

1.005 

Ca-5C|| 

0.54 

0.65 

0.967 

0.971 

Ca -BCl 

0.38 

0.69 

0.979 

0.964 

Ca-AZ?0|| 

-0.29 

-0.33 

0.995 

1.002 

Ca-AZ?0± 

-0  49 

-0  72 

0  993 

1  029 
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Table  II.  Calculated  hydrogen  stretching  modes  in  cm-1  for  the  BCl{  and 
AB01  sites  without  and  with  Ca  and  compared  to  those  in  experiments.  For 
the  calculations,  the  harmonic  frequencies  and  the  anharmonic  contributions 
are  also  shown.  wER  (117  cm-1)  is  the  systematic  error  of  the  calculated 
O-H  vibration  frequencies  (see  text  for  detail). 


Site 

Present 

Measured  (o 

o>0 

Ao> 

w=w0+Ao)+wer 

BC\\ 

3618 

-314 

3421 

3611a 

ABq  ± 

3318 

-326 

3109 

Ca-5C|| 

4022 

-205 

3934 

Ca-AZ?0± 

3443 

-353 

3207 

3326b 

aReference  12. 
bReference  6. 


conductor  with  O  being  the  negative  charge  center  and  Zn 
being  the  positive  charge  center.  Thus,  qualitatively  speak¬ 
ing,  Ca  has  a  stronger  tendency  to  transfer  electron  to  its 
neighboring  O  atoms  than  Zn.  In  other  words,  O  atoms 
neighboring  a  Ca  will  be  more  negatively  charged  than  any 
other  host  O  atoms.  A  more  negatively  charged  O  ion  will 
bind  the  positively  charged  H+  ion  more  strongly.  This  makes 
the  Ca-A^oi  site  most  stable.  Conversely,  a  more  positively 
charged  Ca  ion  will  repel  more  of  the  positively  charged  H+ 
than  any  host  Zn  ion  does.  This  makes  the  Csl-BC\\  sites  least 
stable. 

To  further  our  qualitative  discussion  above,  we  employ  a 
Bader  analysis.22,23  Our  calculations  show  that  in  bulk  ZnO, 
each  Zn  atom  donates  1.13  electrons  to  the  neighboring  O 
atoms.  In  other  words,  Zn  can  be  considered  as  a  positive 
center  with  q=  + 1.1 3e,  whereas  O  can  be  considered  as  a 
negative  center  with  q=- I.13e.  This  explains  why  an  H+ 
prefers  to  stay  close  to  the  O  instead  of  the  Zn.  Bader  analy¬ 
sis  on  CaZn  further  shows  that  Ca  donates  lA2e,  which  is 
0.29e  more  than  that  of  Zn.  The  extra  donation  is  almost 
equally  distributed  among  the  four  O  nearest  neighbors: 
0.06p  to  the  O  atom  in  the  [0001]  direction  and  0.07 p  to  each 
of  the  remaining  three  O  atoms  in  the  perpendicular  direc¬ 
tions.  Thus,  charge  transfer  from  Ca  to  its  neighboring  O 
atoms  is  clearly  demonstrated. 

Next,  we  need  to  determine  with  the  method  in  Ref.  13  if 
the  new  configurations  would  affect  the  O-H  stretching  fre¬ 
quencies  and  the  results  are  shown  in  Table  II.  Although  the 
absolute  values  here  differ  from  experiments  by 
100-200  cm-1,  we  find  that  the  difference  between  the  two 
experimental  modes  (3611  cm_1-3326  cm-1  =  285  cm-1)  is 
in  reasonable  agreement  with  that  between  the  Ca-free  BC\\ 
site  and  Csl-ABq1  site  in  the  calculation 
(3421  cm-1-3207  cm-1  =  214  cm-1).  Our  calculated 
Ca-AB01  oscillator  is  oriented  107°  away  from  the  c-axis, 
which  is  consistent  with  the  experiment  on  the  Cermet 


samples  (with  a  112°  angle),  whereas  our  calculated  BC\\ 
oscillator  is  along  the  c-axis,  which  is  also  consistent  with 
experiments  on  the  EaglePicher  samples. 

IV.  CONCLUSION 

First-principles  calculations  reveal  that  the  isovalent  CaZn 
in  ZnO  has  an  unexpected  strong  Coulomb  binding  with  H+ 
ion  on  the  Ca -AB01  site.  Our  finding  is  consistent  with  re¬ 
cent  experiments  and  may  explain  the  controversies  on  the 
measurements  on  the  Cermet  samples.  The  unexpected  bind¬ 
ing  between  Ca  and  H+  is  explained  by  a  new  charge  transfer 
mechanism,  which  could  be  universal  to  ionic  solids  but  ab¬ 
sent  in  covalent  materials. 
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